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random  variable  defined  as  tl^  sun  of  independent  H-f unction 
variates.  This  allows  an  analyst  to  find  the  distribution  of 
more  conplicated  algebraic  combinations  of  independent  random 
variables.  The  method  and  inplementing  conputer  program  sure 
demonstrated  through  five  exanples.  For  conparison,  the  exact 
distribution  of  the  general  sun  of  independent  Erlang  variates 
with  different  scale  parameters  is  derived  using  Laplace 
transforms  and  partial  fractions  deconposition. 

The  H-function  is  the  most  general  special  function, 
encompassing  as  a  special  case  nearly  every  named  mathematical 
function  and  continuous  statistical  distribution.  The  Laplace 
and  Fourier  transforms  (and  their  inverses)  and  the  derivatives 
of  an  H-function  are  readily-determined  H-functions.  The 
Mellin  transform  of  an  H-function  is  also  easily  obtedned.  The 
H-function  exactly  represents  the  probability  density  fimction 
and  cuiulative  distribution  function  of  nearly  all  continuous 
statistical  distributions  defined  over  positive  values. 

A  previously  unstated  restriction  on  the  variable  in  the 
H-function  representations  of  power  functions  and  beta-type 
functions  is  highlighted.  Several  ways  of  overcoming  this 
limitation  vdien  representing  mathematical  functions  are 
presented.  The  restriction,  however,  is  an  advantage  when 
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r^resenting  certain  statistical  distributions.  Many  new  H- 
function  representations  of  other  nethenatical  functions  are 
also  given. 

The  hierarchical  structure  among  classes  of  H-f unctions  is 
given  through  seven  new  theorens.  Every  class  of  H-f  unctions 
is  wholly  contained  in  nsny  higher-order  classes  of  H-functioois 
throiigh  the  application  of  the  diq>lication,  triplication,  and 
nultiplication  fomulas  for  the  ganna  function. 

Four  new  theorems  show  when  and  how  a  generalizing 
constant  may  be  present  in  an  H-f  motion  representation.  Many 
generalized  H-fmction  representations  are  given,  including 
those  of  every  cumilative  distribution  f motion  of  an  H- 
f  met  ion  variate. 
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CHAPTER  1 


INTROmCTION  AND  REVIEM 

1.1.  PURPOSE  AND  SCOPE 

The  primary  purpose  of  this  research  effort  was  to  develop 
a  practical  method  of  finding  an  H-functicoi  distribution  for 
the  sun  of  two  or  more  independent  H- function  variates.  Sinple 
formulas  exist  vdiich  imnediately  give  the  probability  density 
function  (p.d.f.),  as  an  H-f unction  distribution,  of  the  random 
variable  defined  as  the  product,  quotient,  or  power  of 
independent  H-f  unction  variates.  Unfortunately,  there  are  no 
similar  formulas  for  the  sxm  or  difference  of  independent  H- 
f unction  variates. 

A  related  issue  was  whether  the  class  of  H-functions  is 
closed  under  the  operation  of  multiplication.  In  other  words, 
is  the  product  of  two  H-functions  another  H-function?  It  is 
iirportant  to  make  the  distinction  here  between  the  product  of 
two  H-functions  and  the  p.d.f.  of  the  random  variable  defined 
as  the  product  of  two  H-function  variates.  It  is  well  known 
that  the  latter  case  is  an  H-function.  But  the  former  case  was 
unproven.  Of  course,  similar  statements  can  be  made  about  the 
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quotient  of  two  H-f unctions. 

If  the  class  of  H-functions  is  closed  vmder 
nul tipi i cation,  one  could  easily  find  the  p.d.f.  (as  an  H- 
fvnction)  of  the  svin  or  difference  of  independent  H-fmction 
variates.  The  Laplace  (or  Fourier)  transforms  of  the  H- 
f unction  variates  in  the  sixn  (or  difference)  are  imnediately 
available  as  H-f  unctions  of  higher  order.  The  product  of  these 
H-functions  (in  transform  space)  would  yield  the  transform  of 
the  desired  density.  If  this  product  was  available  as  another 
H-function,  it  could  be  inverted  from  transfonn  space 
analytically. 

Because  the  H-function  can  exactly  represent  nearly  every 
comnan  mathematical  function  and  statistical  density,  there  was 
anple  reason  to  suspect  that  the  product  of  two  H-f\]nctions 
was,  in  general,  another  H-fxmction.  Indeed,  there  are  many 
cases  where  two  individual  fiaicticons  and  their  product  are  all 
special  cases  of  the  H-function. 

Throughout  this  thesis,  a  nunnber  of  other  new  results  are 
identified  with  an  asterisk.  Sufficient  convergence  conditions 
for  the  alternate  definition  of  the  H-function  are  given  in 
Section  2.3.  These  show  how  the  H-function  may  be  evaluated  by 
the  sum  of  residues,  without  first  changing  the  form  of  the 
alternate  definition  of  the  H-function  to  that  of  the  primary 
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definition. 

The  hierarchical  stmctiire  among  classes  of  H-fvmctians  is 
given  through  seven  new  theorems  in  Section  2.4.6.  Every  class 
of  H-f unctions  is  wholly  contained  in  nany  higher-order  classes 
of  H-fmcticns  through  the  application  of  the  di;Q>lication, 
triplication,  and  nultiplication  formulas  for  the  gamna 
function.  Figure  1  in  Section  3.5  illtistrates  this 
hierarchical  structvure  with  a  venn  diagram  showing  many  coninon 
statistical  distributions  as  first  and  second  order  H-function 
distributions. 

Four  new  theorems  in  Section  2.4.7  show  when  and  how  a 
generalizing  constant  nay  be  present  in  an  H-function 
representation.  Many  generalized  H-function  representations 
are  given,  including  those  of  every  cvnulative  distribution 
function  of  an  H-function  variate.  The  generalizing  constant 
is  also  possible  in  the  H-fvnction  representations  of  power 
functions,  the  error  function  and  its  ccnplement,  the 
incomplete  ganna  function  and  its  complement,  the  incomplete 
beta  function  and  its  complement,  many  inverse  trigonometric 
and  hyperbolic  functions,  and  the  logarithmic  functions. 

A  nxarber  of  new  H-fvnction  representations  of  certain 
mathematical  functions  and  statistical  distributions  are  given 
in  Sections  2.5,  2.6,  and  3.6.  Several  of  these  expand  upon  a 
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previously  unstated  limitation  cn  the  variable  in  the  H- 
function  representations  of  power  functions  and  beta-type 
functions . 

The  exact  distribution  of  the  general  sun  of  independent 
Erlang  variates  with  different  scale  parameters,  X,  is  derived 
in  Section  5.1.3.  An  Erlang  variate  is  simply  a  gamma  variate 
with  an  integer  shape  parameter  r.  The  derivation  uses  partial 
fractions  to  decompose  the  product  of  Laplace  transforms  of  the 
individual  densities.  This  produces  a  sun  of  terms,  each  of 
which  can  easily  be  inverted  front  transform  space,  yielding  the 
desired  density  of  the  sun  of  independent  Vcuriates. 

Since  the  H-fimction  is  not  defined  for  zero  or  negative 
real  argunents,  the  scope  of  this  research  effort  was  limited 
to  continuous  random  variables  defined  only  over  positive 
values.  Continuous  and  doubly  infinite  distributions  such  ats 
the  normal  and  Student's  t  are  only  represented  as  H-functions 
in  their  folded  forms. 

1.2.  LITERATURE  SURVEg 

Regrettably,  little  research  in  the  field  of  H-functioons 
has  been  done  in  the  Uhited  States.  Much  of  what  is  known 
about  the  H-f unction  is  due  to  Indian  mathematicians.  Mathai 
and  Saxena  [1978]  and  Srivastava  et  al  [1982]  compiled  many 
results  of  the  early  study  of  H-fxmctions.  In  recent  years. 
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Soviet  nathematicians  [Prxidnikov  et  al ,  1990]  have  shown  a 
considerable  interest  in  the  H-fvmction  and  have  developed 
significant  new  results. 

The  foundation  of  H-ftmction  theory  is  grounded  in  the 
gamna  function,  integral  transform  theory,  conplex  analysis, 
and  statistical  distribution  theory.  Therefore,  several 
landmark  references  such  as  Abramowitz  and  Stegxai  [1970], 
Erd^lyi  [1953],  Erdelyi  [1954],  and  Springer  [1979],  though 
scmetdiat  dated,  have  timeless  value. 

Carter  [1972]  defined  the  H-function  distribution  and, 
using  Mel  1  in  transform  theory,  gave  startling  and  powerful 
results  showing  that  products,  quotients,  and  rational  powers 
of  independent  H-function  variates  are  themselves  H-function 
variates.  Further,  the  p.d.f.  of  the  new  random  variable  can 
imnediately  be  written  as  an  H-function  distribution.  Ilie 
visual  techniques  of  conditioning  on  one  of  the  random  variables 
and/or  using  the  Jacobian  of  the  transformation  are  no  longer 
necessary. 

The  above  results  become  especially  useful  vdien  one 
realizes  that  nearly  every  comnnon  positive  continuous  random 
variable  can  be  written  as  an  H-function  distribution. 
Therefore,  the  p.d.f.  of  any  algebraic  combination  involving 
products,  quotients,  or  powers  of  any  nvnber  of  independent 


positive  continuous  random  variables  can  imnediately  be  written 
as  an  H-function  distribution. 

Carter  [1972]  also  wrote  a  conrputer  program  to 
calculate  the  moments  of  an  algebraic  combination  of 
independent  H-function  variates  and  approximate  the  p.d.f.  and 
cumulative  distribution  function  (c.d.f.)  from  these  moments. 
The  approximation  procedure  was  developed  by  Hill  [1969]  and, 
if  possible,  uses  either  a  Gram-Charlier  type  A  series  (Hermite 
polynomial)  or  a  Laguerre  polynomial  series.  If  a  series 
approximation  is  not  possible,  the  first  four  moments  are  used 
to  fit  a  probability  distribution  from  the  Pearson  family.  As 
Carter  [1972]  himself  notes  "...  there  were  many  situations  in 
which  the  methods  did  not  work  or  in  idiich  the  approximations 
were  totally  unsatisfactory." 

Springer  [1979]  literally  wrote  the  book  on  the  algebra  of 
random  variables.  He  gives  an  excellent  explanation  of  the 
value  of  integral  transforms  in  finding  the  distribution  of 
algebraic  combinations  of  random  variables.  He  also  gave  the 
known  applications  of  the  H-fmction  in  these  problems. 

Cook  [1981]  gave  a  very  thorough  survey  and  an  extensive 
bibliography  of  the  literature  related  to  H-functions  and  H- 
function  distributions.  He  also  presented  a  technique  for 
finding,  in  tabular  form,  the  p.d.f.  and  c.d.f.  of  an  algebraic 
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combination  (including  sixns  and  excluding  differences)  of 
independent  H-f unction  variates. 

Cook's  technique  [1981;  Cook  and  Barnes,  1981]  first  uses 
Carter's  [1972]  results  to  find  the  H-f unction  distribution  of 
any  products,  quotients,  or  powers  of  random  variables.  The 
Laplace  transform  of  each  term  in  the  resulting  sum  of 
independent  H-f unction  variates  is  then  obtained.  These 
transform  functions  are  evaluated  and  multiplied  at 
corresponding  values  of  the  transform  variable,  yielding  a 
tabular  representation  of  the  Laplace  transform  in  transform 
space.  This  Laplace  transform  is  then  nixnerically  inverted 
using  Cnitp's  method.  His  PCSiTRAN  computer  program  implements 
this  technique  and  will  plot  the  resulting  p.d.f .  and  c.d.f . 

Bodenschatz  and  Boedigheimer  [1983;  Boedigheimer  et  al, 
1984]  developed  a  technique  to  fit  the  H-f  unction  to  a  set  of 
data  using  the  method  of  moments.  The  technique  can  be  used  to 
curve-fit  a  mathematical  function  or  to  estimate  the  density  of 
a  particular  probability  distribution.  Their  FCSTTRAN  computer 
program  will  accept  known  moments,  univariate  data,  ordered 
pair  data  from  a  relative  frequency,  or  ordered  pair  data 
directly  from  the  function. 

Kellogg  [1984;  Kellogg  and  Barnes,  1987;  Kellogg  and 
Barnes,  1989]  sttxiied  the  distribution  of  products,  quotients. 
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axid  powers  of  dependent  random  variables  with  bivariate  H- 
function  distributions.  Jacobs  [1986;  Jacobs  et  al,  1987] 
presented  a  method  of  obtaining  parameter  estimates  for  the  H- 
functicn  distribution  using  the  method  of  naxiirun  likelihood  or 
the  method  of  moments. 

Prudnikov  et  al  [1990]  gave  extensive  tables  of  H-f unction 
results  and  Mellin  transforms.  Their  books,  although  more 
terse  than  the  series  by  Erdelyi  [1953  and  1954],  are  at  least 
as  complete  and,  likely,  will  become  the  new  standard  reference 
for  special  functions. 

1.3.  INTE}GRRL  TRftNSFCiRMS  AND  TRANSFORM  PAIRS 

Integral  transforms  are  frequently  encountered  in  several 
areas  of  mathematics,  probability,  and  statistics.  Although 
various  integral  transforms  exist,  certain  characteristics  are 
corrmon  among  them.  The  function  to  be  transformed  is  usually 
multiplied  by  another  function  (called  the  kernel)  and  then 
integrated  over  an  appropriate  range.  What  distinguishes  the 
various  transforms  are  the  kernel  fmction,  the  limits  of 
integration,  and  the  type  of  integration  (e.g.  Riemann  or 
Lebesgue) . 

Often  the  use  of  integral  transforms  can  simplify  a 
difficult  problem.  Laplace  transforms  are  usually  first 
encountered  in  the  solution  of  systems  of  linear  differential 
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equations.  In  probability  and  statistics,  certain  integral 
transfonns  are  often  known  by  other  names  such  as  the  mcment 
generating  function,  characteristic  function,  or  probability 
generating  f motion. 

The  definitions  of  mcjst  integral  transforms  are  not 
standard.  It  is  important,  therefore,  to  explicitly  state  the 
form  of  the  definition  to  be  used.  Listed  below  are  the 
definitions  of  certain  integral  transforms  and  the 
corresponding  inverse  transforms  as  used  in  this  thesis. 
Together,  each  transform  and  its  corresponding  inverse 
constitute  a  transform  pair. 

1.3.1.  LftPIACE  TRANSFORM 

Consider  a  f motion  f(t)  which  is  sectional ly  continuom 
and  defined  for  all  positive  values  of  the  variable  t  with 
f(t)=0  for  tiO.  A  sectional ly  continuous  f motion  may  not  have 
an  infinite  nixrber  of  discontinuities  nor  any  positive  vertical 
asymptotes.  If  f(t)  grows  no  faster  than  an  exponential 
fmction,  then  the  Laplace  transform  of  f(t)  will  exist.  There 
must  exist  two  positive  numbers  M  and  T  such  that  for  all  t>T 
and  for  some  real  number  a. 


(1.1) 
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The  definition  of  the  Laplace  transform  of  the  function 


f(t),  f(t)  is 

*s{  '  f 

J  n 


e  f(t)  dt 


(1.2) 


In  general,  s  is  a  ccnplex  variable.  The  Laplace  transform  of 
f(t)  vri.ll  exist  for  the  real  part  of  s  greater  than  a 
(Re(s)  >  a). 

The  inversion  integral  or  inverse  Laplace  transform  is 
given  by 

w+i® 

I  *s{ 

'*w-ioo 

where  f(t)  ^  is  an  analytic  function  for  Re(s)  >  w.  A 

function  is  analytic  at  s=s  if  its  derivative  exists  at  s  and 

o  o 

at  every  point  in  some  neigMx>rhood  of  s^.  The  Taylor  series 
expansion  of  an  analytic  function  of  a  ccnplex  variable  will 
exist,  converge,  and  equal  the  function  evaluated  at  the 
arg\xnent.  For  all  practical  purposes,  a  function  f(t)  and  its 
Laplace  transform  (if  it  exists)  uniquely  determine  each  other. 

In  probability  and  statistics,  if  f(t)  is  the  p.d.f.  of  a 
random  variable  defined  only  for  positive  values,  its  moment 
generating  fixnction  is  siirply  the  Laplace  transform  with  r 
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replacing  -s  in  Eg  (1.2). 

1.3.2.  FOURIER  'HaNSFORM 

The  form  of  the  exponential  Fourier  transform  used  in  this 
thesis  is 


f(t)  J  =  e^®^  f(t)  dt  (1.4) 

~a> 

The  Fourier  transform  is  a  fmction  of  the  ccnplex  variable  s. 

If  f(t)  is  a  p.d.f.,  this  definition  corresponds  to  the 
definition  of  the  characteristic  function  in  probability  and 
statistics.  The  characteristic  function  of  a  p.d.f.  will 
always  exist  but  the  mcment  generating  function  of  a  p.d.f.  may 
or  may  not  exist. 

The  inversion  integral  or  inverse  Fourier  transform  is 
given  by 


vdiere 


lim  f(t)  +  lim  f(t) 

tr*t  fc+t^ 

o  o 

t>t„ 

o  o 


(1.6) 


If  f(t)  exists  and  is  continxious  at  t  ,  the  inverse  Fourier 

o 
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transform  of  f(t)  j  v*ill  give  f (tj  •  If  f(t)  is  not 

continvious  at  t^,  the  inverse  Fourier  transform  of  5^^  f(t)  ^ 

will  produce  the  average  of  the  limits  of  f(t)  fron  the  left  of 

t  and  the  right  of  t  . 
o  o 


1.3.3.  MELLIN  TRANSFCMi 

Becaiise  the  Mel  1  in  transfonn  is  perhaps  less  well  known 
than  the  Laplace  or  Fourier  transforms  and  because  the  Mellin 
transform  is  so  crucial  in  the  study  of  H-f unctions,  both 
ccrnnon  sets  of  transform  pairs  will  be  presented.  The  Mellin 
transform  uses  a  power  function  instead  of  an  exponential 
function  as  its  kernel. 

Again  consider  a  fimction  f(t)  which  is  sectionally 
continuous  and  defined  for  all  positive  values  of  the  variable 
t  with  f(t)=0  for  tsO.  Using  what  will  be  regarded  in  this 


thesis  as  the  primary  definition  of  the  Mellin  transform,  the 
Mellin  transform  of  f(t),  it  /  f(t)  \,  is 


^s{ 

'*0 


(1.7) 


The  Mellin  transform  is  related  to  the  Fourier  and  Laplace 
transforms  as  follows  [Erd^lyi,  1954,  p.305]: 
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M  *(•*)  }  (1-8) 

'  *-s{  ‘C®'*)  } 

Again,  s  is  a  conplex  variable.  The  Mellin  transform 
inversion  integral,  or  inverse  Mellin  transform,  is  given  as 
(i)+iao 

~  ^  I  ^s{  ^  (1-10) 

“^w-ioo 

As  mentioned  earlier,  there  is  another  inportant  transform 
pair  also  referred  to  as  a  Mellin  transform  pair.  This 
alternate  definition  will  arise  later  tdmn  an  alternate 
definition  of  the  H-functicn  is  given.  The  alternate 
definition  is 

P® 

it'l  f(t)  J  f(t)  dt  (1.11) 

Jq 

with  inverse  transform 

(t»+ia> 

il'l  f(t)  I  dr  (1.12) 

p-i® 

1.4.  TRANSFQBMATIOHS  OF  INDEPEMEOTT  RANDOM  VARTART.RS 


A  conmon  problem  in  statistical  distribution  theory  is  to 
find  the  distribution  of  an  algebraic  conbination  of 
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indQ>endent  randan  variables.  The  algebraic  conbination  could 
include  suns,  differences,  products,  and/or  quotients  of 
independent  randan  variables  or  their  powers.  It  is  iirportant 
to  recognize  that  the  algebraic  conbination  is  itself  a  random 
variable  and,  therefore,  has  a  probability  distribution.  The 
task  is  to  find  this  distribution. 

Using  the  properties  of  nathematical  expectation,  it  is 
relatively  easy  to  find  the  mean,  variance,  and  other  monents 
of  the  algebraic  corbination  of  independent  random  variables. 
For  example,  the  mean  of  the  sun  of  two  indepeafvient  razxicm 
variables  is  simply  the  sun  of  the  means  of  the  random 
variables.  Finding  the  completely  specified  distribution  of 
the  algebraic  conbination  is  usually  much  more  difficult. 

For  simple  combinations  of  independent  random  variables, 
the  method  of  Jacobians  is  often  employed.  An  appropriate 
one-to-one  transformation  between  the  independent  random 
variables  in  the  algebraic  conbination  and  a  set  of  new  random 
variables  is  first  created.  After  finding  the  inverse 
transformation  fvnctions,  the  Jacobian  (the  determinant  of  the 
matrix  of  first  partial  derivatives  of  the  inverse  functions) 
may  be  computed. 

The  joint  p.d.f .  of  the  newly  defined  random  variables  is 
the  absolute  value  of  the  Jacobian  multiplied  by  the  product  of 
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the  original  densities  with  the  inverse  transfornaticn 
f motions  substituted  for  the  variables.  A  great  deal  of  care 
nust  be  used  in  determining  the  values  of  the  new  variables  for 
which  the  joint  density  is  nonzero.  Once  this  is  done,  the 
desired  marginal  density  can  be  obtained  by  integrating  over 
the  carpi ete  ranges  of  the  other  new  variables. 

An  exanple  of  the  method  of  Jacobians  will  be  presented 
below.  In  most  of  the  following  sections,  however,  only  the 
result  using  integral  transforms  will  be  given. 


1.4.1.  DISTRIBOTIOM  CF  A  SOM 

Let  and  Xj  be  independent  random  variables  with 
respective  densities  and  f2(*2) '  nonzero  only  for 
positive  values  of  the  variable.  Svppose  we  want  the  density 
of  y=Xj^+X2.  Using  the  method  of  Jacobians,  we  define  WsXj  so 
the  inverse  transformations  are  Xj^=Y-W  and  X2=W.  The  Jacobian 
is 


J  = 


1 

0 


(1.13) 


The  joint  density  of  Y  and  W  is 

f^(y,w)  =  fj^(y-w)  fjCw)  0  <  w  <  y  <  ®  (1.14) 

The  marginal  p.d.f.  of  Y  can  be  obtained  by  integrating  the 
joint  p.d.f.  with  respect  to  H  over  the  range  of  W. 
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t^iy)  =  J  ^ 


(1.15) 


Eq  (1.15)  inay  be  recognized  as  the  Fourier  convolution 
integral  [Springer,  1979,  p.  47].  This  is  no  accident  or 
coincidence  as  the  Fourier  or  Laplace  transform  could  also  have 
been  used  to  find  the  distribution  of  Y.  It  is  well  known  that 
the  Laplace  (or  Fourier)  transform  of  the  p.d.f .  of  the  stn  of 
two  independent  random  variables  defined  for  positive  values  is 
the  product  of  the  Laplace  (or  Fourier)  transforms  of  the 
individual  densities.  Further,  the  product  of  two  transfonn 
functions,  upm  inversion,  yields  a  convolution  integral. 

If  the  product  of  transfonn  functions  can  be  recognized  as 
the  transform  of  some  function,  then  the  convolution  inversion 
is  not  necessary.  The  p.d.f.  of  Y  is  the  fmction  >diose 
transform  is  the  product.  Nhen  statisticians  use  the  moment 
generating  function  of  each  density  to  find  the  p.d.f.  of  Y, 
this  recognition  approach  is  usually  taken.  One  advantage  of 
vising  transform  functions  is  that  the  procedure  easily  extends 
when  the  distribution  of  the  sun  of  three  or  more  independent 
random  variables  is  desired. 

Finding  the  distribution  of  the  sun  of  independent  random 
variables  with  certain  special  distributional  forms  is 
considerably  sinplified.  Several  of  these  cases  are  covered  in 
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the  following  stibsections . 


1.4. 1.1.  INFINITE  DIVISIBILITY 
Although  a  connplete  discussion  of  infinite  divisibility  is 
beyond  the  scope  of  this  thesis,  its  definition  is  given  below 
[Petrov,  1975,  p.  25]. 

A  distribution  fxmction  F(x)  and  the  corresponding 

characteristic  function  f(t)  are  said  to  be 

infinitely  divisible  if  for  every  positive  integer  n 

there  exists  a  characteristic  function  f  (t)  such 

n 

that 

f(t)  =  (yt))" 

In  other  words,  the  distribution  F  is  infinitely 
divisible  if  for  every  positive  integer  n  there 

★n 

exists  a  distribution  function  F^  such  that  F=F^  . 

n  n 

Here  F**'  is  the  n-fold  convolution  of  the  function 


Coninan  exanples  of  infinitely  divisible  distributions  include 
the  nomal  and  Poisson  distributions. 


1.4. 1.2.  SPECIAL  CASES 

The  distribution  of  the  stin  of  independent  random 
variables  with  certain  distributional  forms  is  well  known  and 
iirmediately  available.  For  example,  the  sunn  of  n  independent 
and  identically  distributed  random  variables  with  a  Bernoulli 
distribution  with  parameter  p  has  a  Binomial  distribution  with 
parameters  n  and  p.  Similarly,  the  sun  of  independent 
geometrically  distributed  random  varied)!  es  with  a  common 


18 


parameter  has  a  negative  Bincndal  (or  Pascal)  distribution. 

In  continuous  random  variables,  the  sixn  of  n  independent 
and  identically  distributed  random  variables  vrith  an 
exponential  distribution  with  parameter  X  has  a  ganma 
distribution  with  parameters  n  and  X. 

1.4. 1.3.  REPRODUCTIVE  mffTRlBIJTlOHS 
A  probability  distribution  is  "reproductive"  if  it 
replicates  under  positive  addition  of  independent  random 
variables  with  the  same  distributional  form.  The  noxnBl 
distribution  is  reproductive  since  given  that 

-  Normal  ^  *^2'  ‘^l  )  ^  ~  Normal  ^  M2»  <^2  ) 

and  X^  are  statistically  independent,  and 

Y  =  Xj^  +  Xj 

then 

Y  ~  Normal  ^  ] 

Other  examples  of  reproductive  distributions  include  the  Chi- 
Square  and  Poisson  distributions. 

It  is  well  known  that  the  ganma  distribution  is 
reproductive  provided  the  scale  parameter,  X,  is  the  same  for 
each  random  variable  in  the  sutn.  In  particular,  if 

-  Ganma  [  *^2'  ^1  ]  ^  ~  Ganma  f  r2,  Xj  1 


and  X2  are  statistically  independent, 
and  y  =  +  Xj 
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then 


Y  ~  Ganma 


(  ‘  ) 


The  result  for  ganma  distributions  above  is  readily 
verified  by  considering  the  product  of  the  Laplace  transforms 
of  each  p.d.f.  Using  the  definition  in  Eg  (1.2),  the  Laplace 
transform  of  the  ganma  p.d.f.  with  parameters  r  and  X  is 

[  "3^"]  Clearly,  then,  if  X^  and  X^  are  independent  random 

variables  with  ganma  distributions  and  a  conmon  scale 
parameter,  X, 


(1.16) 

(1.17) 

f  X  1V‘^2 
-  l-s+rj 

(1.18) 

which  is  recognized  as  the  Laplace  trazvsform  of  a  ganma  p.d.f. 
with  parameters  confirms  the  reproductive 
property  for  the  ganma  distribution  when  X  is  conmon. 


1.4.2.  DISroiBUTIOH  OF  A  DIFFBREMCE 

The  study  of  the  distribution  of  the  difference  of 


independent  random  variables  has  received  nuch  less  attention 
than  for  the  sun.  The  normal  distribution  is  one  notable 


exception  since  given  that 

-  Normal  ^  j  ^  "  Normal  ^  /ij, 

and  are  statistically  independent,  and 
Y  =  X^  -  X2 

then 


Y  ~  Normal  ^  ] 

If  X^  and  X2  are  independent  random  variables  with 
respective  densities  f 3“^  ^2(^2)  '  nonzero  only  for 
positive  values  of  the  variable,  then  the  density  of  YsX^^-Xj  is 
the  inverse  Fourier  transform  [^ringer,  1979,  p.  59] 


for  -w  <  y  <  a> 


ds 


(1.19) 


1.4.3.  DISTOIBUTIOW  OF  A  PROTOCT 

If  X^  and  X2  are  iiulependent  random  variables  with 
respective  densities  and  £2(^2)  '  nonzero  only  for 

positive  valxies  of  the  variable,  then  the  density  of  Y=(X^)  (X2) 
is  the  inverse  Mellin  transform  [Springer,  1979,  p.  97] 
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(i)+100 

^  ^  y  ^  ^s{  *lC*l)  }  ^s{  *2  W  }  ^ 

'*w-ioo 


for  0  <  y  <  00 


(1.20) 


1.4.4.  DISTRIBUTIOM  OF  A  OCIOTIEMT 

If  and  are  independent  random  variables  with 
respective  densities  ^2(^2)  '  nonzero  only  for 

positive  values  of  the  variable,  then  the  density  of  —  is 

the  inverse  Mel 1 in  transform  [Springer,  1979,  p.  100] 
u+ioo 

I  *3{  }  *2-s{  ‘jW  }  ^ 

''w-ioD 

for  0  <  y  <  00  (1.21) 

One  conmon  example  of  the  quotient  of  two  independent 
random  variables  lies  in  tl%  derivation  of  the  Snedecor  F 
distribution  [Springer,  1979,  pp.  328-9].  If 

Xj^  -  Chi-Square  iu)  ~  Chi-Square  (u) 

X^  and  X2  are  statistically  independent,  and 


Y  ~  P  (p,«) 
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1.4.5.  DISTRIBOTICTl  OF  A  VaRIATE  TO  A  POWER 
If  X  is  a  continiious  random  varicdble  with  density  £(x), 
nonzero  only  for  positive  values  of  the  variable,  then  the 
density  of  Y=X^  is  the  inverse  Mel  1  in  transform  [Springer, 
1979,  p.  212] 

tj+ioo 

''o-iflo 

for  0  <  y  <  00  (1.22) 

Gne  comnnon  exarrple  of  a  variate  to  a  power  is  finding  the 

distribution  of  the  square  of  a  standard  normal  (zero  mean, 

unit  variance)  random  variable.  If  X  -  Normal  (0,1)  and  Y=X^ 

then  Y  ~  Chi-Square  (1).  The  same  result  holds  if  X  has  a 

2 

half-normal  distribution  with  a  =1  [Springer,  1979,  pp.  213-4]. 


1.4.6.  MOMENTS  OF  A  DISTRIBUTIOW 

If  X  is  a  continuous  random  variable  with  density  f(x), 
nonzero  only  for  positive  values  of  the  variable,  then  the 
moments  about  the  origin  of  f(x)  are 

=  e[  ] 


(1.23) 


provided  the  integral  in  Eq  (1.23)  exists. 
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There  is  a  natural  relationship  between  the  integral  given 
in  Eq  (1.23)  and  the  Mellin  transform  of  the  function  f(x) 
given  in  Eq  (1.7)  or  Eq  (1.11).  Using  the  Mellin  transform  in 
Eq  (1.7)  we  can  write 

^r  ’  (1-24) 

This  relationship  sirnplifies  the  coninitation  of  mcments  for  H- 
f motions  and  H-f motion  distributions. 


CHAFTEK  2 


TOE  H-FTJNCTIC3N 


The  H-fmction  is  a  very  general  function,  enccnpassing  as 
special  cases  nearly  every  named  mathenatical  function  and 
continuous  statistical  distribution  defined  over  positive 
values.  Although  the  H-f unction  does  not  enjoy  an  extensive 
popularity  and  acceptance  in  the  fields  of  nathenetics, 
probability,  and  statistics,  this  is  primarily  becatise 
mathematicians  and  statisticians  have  not  yet  learned  of  its 
versatility  and  power.  Most  analysts  are  not  familiar  with  the 
H-f  unction,  and  many  who  have  seen  the  H-f  unction  definition 
may  have  been  disquieted  by  its  overt  complexity.  TOis  is 
unfortunate  because  practical  use  of  the  H-function  does  not 
require  extensive  knowledge  of  complex  analysis  and  integral 
transform  theory. 

Mathematical  functions  defined  by  an  integral  tdiich  do  not 
have  a  closed  form  representation  are  comnon.  Examples  include 
the  ganrna  function,  the  error  function,  and  the  cunnulative 
normal  probability  density  fmction.  In  all  of  these  cases, 
the  function  is  usually  evaluated  with  the  help  of  tables.  The 
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H-fmction  is  another  exairple  of  such  functions.  Like  all  of 


the  transcendental  functions  (e.g.  e*,  sin  x,  cos  x),  it  must 
be  evalviated  using  an  infinite  series  expansion.  A  FORTRAN 
ccnputer  program  is  available  [Cook,  1981;  Cook  and  Barnes, 
1981]  which  will  evaluate  the  H-fvinction  at  desired  values  once 
the  parameters  are  specified. 

Fox  [1961]  first  developed  the  H-function  as  a  direct 
generalization  of  Meijer's  G-f unction.  Mathai  and  Saxeaia 
[1978]  presented  many  \iseful  properties  of  the  H-function  and 
listed  the  mathematical  functions  that  are  known  to  be  special 
cases  of  the  H-function.  More  recently,  Prudnikov  et  al  [1990] 
conpiled  e.'ttensive  results  of  all  special  ftmctions,  including 
H-f unctions. 

2.1.  PRIMARy  ICTINITION 

The  primary  definition  of  the  H-function  as  used  in  this 
thesis  is: 


H(«)  =  ”1 


cz 


i=l,...,p  j=l,...,q 


1 

2xi 


C 


m 

n  r 
j=i 

p 

n  r 

i=n+l 


n 

n  r 

i=l 

q 

n  r 
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f^re  z,  c,  and  all  and  are  real  or  conplex  nvirbers,  all 
and  are  positive  real  nizibers,  and  m,  n,  p,  and  q  are 
integers  such  that  O^n&q  and  Qsnsp.  Ehpty  products  are  defined 
to  be  vnity  (1).  The  path  of  integration,  C,  is  a  contour  in 
the  carplex  s-plane  fran  u-ico  to  u-fico  such  that  all  Left  Half- 

Plane  (lilP)  poles  of  n  rfb.-t-B-sl  lie  to  the  left  of  C  and  all 

j=l  1  3  3  J 

Right  Half-Plane  (RHP)  poles  of  n  T^l-a^-A^sj  lie  to  the 
right . 

Although  the  H- function  is  defined  for  the  conplex 

varis^le  z,  vte  viill  often  restrict  our  attention  to  the  real 

variable  x.  Further,  since  the  H-function  is  not  valid  for 

non-positive  real  values  of  z,  we  will  often  consider  only 

positive  values  of  the  real  variable  x. 

The  definition  of  the  H-function  in  Eg  (2.1)  nay  be 

recognized  as  the  inverse  Mellin  transform  where  the  transform 

pair  is  as  given  by  Eq  (1.7)  and  Eq  (1.10). 

H-functions  are  sometimes  classified  according  to  their 

order,  the  nurber  of  ganma  terms  in  the  integrand  (p>q),  aixl 

the  placement  of  those  terms.  References  will  be  made  to 

certedn  classes  of  H-functions  as  with  particular  values 

P  q 

for  the  parameters  m,  n,  p,  and  q. 
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2.2.  ALTERNRTE  DEFINITIOM 

It  should  be  noted  that  there  is  an  alternate,  but 
ecpii valent  definition  of  the  H-f unction: 


H(cz)  =  hJ  “ 
P  <3 


i— 1, ...  ,p  j~l-,  •  •  •  ,*3 


1 

2>ii 


li  I 

n  rfa.-A.rl  n  rfi-b.+B 
i=n+l  ^  ■'  j=nH-l  ^  ^ 


(cz)*^  dr 


(2.2) 


vdiere  the  restrictions  on  the  variable  and  parameters  are  as 
above.  Here,  the  path  of  integration,  C' ,  is  a  contour  in  the 
ccnplex  r-plane  from  p-Icd  to  p+ioo  such  that  all  RHP  poles  of 

n  rfb  .-B  .rl  lie  to  the  ri^t  of  C'  and  all  UIP  poles  of 
j=l  V  J  3  J 

n  rfl-a.+A.rl  lie  to  the  left. 
i=l  I  ^  ^  J 

Eg  (2.2)  above  can  be  derived  from  Eg  (2.1)  hy  the 
substitution  r=-s.  Under  this  transformation,  since  the 

substitution  r=-s  rotates  about  the  imaginary  axis  all  poles  of 
the  integrand  and  the  contour,  C,  into  a  contour,  C' ,  idiicdi 
also  separates  the  new  LHP  poles  front  the  new  RHP  poles  in  the 
standard  manner  and  direction.  Again,  the  definition  in 
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Eq  (2.2)  corresponds  to  an  inverse  Mel  1  in  transform  where  the 
transform  pair  is  as  given  by  Eg  (1.11)  and  Eg  (1.12). 

2.3.  SUFFICIENT  CX^tVERGEMCE  a^TOITiqiS 

Cook  [1981;  Cook  and  Barnes,  1981]  gave  sufficient 
conditions  for  which  the  H-function  could  be  evaluated  as  the 
sum  of  residues  in  the  appropriate  half -plane  for  certain 
values  of  the  variable.  In  an  unpublished  working  paper, 
Eldred  et  al  [1979]  applied  the  well-known  convergence 
restrictions  of  Mellin-Bames  integrals  to  the  H-function. 
These  restrictions  were  originally  developed  by  Dixon  and 
Ferrar  [1936]  and  are  also  given  by  Erdelyi  [1953,  Vol.  1]  and 
Prudnikov  et  al  [1990]. 

These  conditions  help  determine  how  to  evaluate  the  H- 
f unction  in  Eg  (2.1)  as  the  sum  of  residues  and  give  the  values 
of  the  complex  variaJsle  z  for  which  this  evaluation  is  valid. 
These  convergence  conditions  specify  the  restrictions  on  the 
argument,  arg(z),  and  the  moduli;is,  |z|,  of  z  to  guarantee 
convergence.  They  indicate  that  the  sum  of  residues  for  the  H- 
function  is  always  convergent  for  positive  real  values  of  the 
variable  and  sometimes  convergent  for  values  of  z  with  a 
nonzero  imaginary  part. 

For  the  H-function  defined  by  Eg  (2.1),  the  convergence 
conditions  are  based  on  the  values  of  D,  E,  L,  and  R  defined 


(2.3) 


D  = 


y  A.  +  y  B .  -  y  a.  -  v  b  . 
i=i  ^  3=1  ^  i^+i  ^  j=4fi  3 


E  = 


p  q 

i=l  ^  3=1  ^ 


L  =  Re 


^  .  q  . 


T  b.--^  -  T  a.  +  -4- 

Ld^i  ^  ^  i=i  "  ^ 


(2.4) 


(2.5) 


(2.6) 


Once  these  values  are  determined,  the  H-£unction  in  Eq  (2.1) 
nay  be  evaluated  by  the  positive  sum  of  ISP  residues,  the 
negative  sum  of  RHP  residues,  or  both,  depending  cn  the  value 
of  the  variable  z.  These  criteria  are  based  ypcn  vdiich  said- 
circle  satisfies  the  hypotheses  of  Jordan's  Lemta.  The  types 
of  convergence  and  the  aq^licable  ranges  for  the  conplex 
variable  z  were  given  by  Eldred  et  al  [1979],  Cook  [1981],  and 
Cook  and  Barnes  [1981]  and  are  repeated  below  in  Table  1: 
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definition  (2.2).  By  rewriting  Eg  (2.2)  as 


1 

2ni 


”,  r(l-a,.A,r)  r(b  .-B  .r) 

n  rfl-b  +B  r]  n  rfa  -A.r] 
j=mH  ^  J  J  J  i=n+l  ^  ^  J 


f— 1 

I  C2  J 


and  using  Eq  (2.1)  by  ccaisidering  the  "A”  terms  as  " 
and  the  "B"  terms  as  "A"  terms,  we  have 


Hnml 

q  p 


j=l,...,q  i=l,...,p  J 


Evaluating  D' ,  E' ,  L' ,  and  R'  as  in  Eq  (2.3) 
Eq  (2.6)  yields 

m  n  q  p 

*  D'  =  I  B .  +  y  A.  -  y  B .  -  y  A. 
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(2.7) 
i"  terms 

(2.8) 

through 

(2.9) 

(2.10) 

(2.11) 
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* 


(2.12) 


Ccnparing  Eg  (2.9)  through  Eg  (2.12)  to  Eg  (2.3)  through 


Eg  (2.6), 

D'  in  Eg  (2.9)  =  D  in  Eg  (2.3)  (2.13) 
E'  in  Eg  (2.10)  =  (-1)  [  E  in  Eg  (2.4)  ]  (2.14) 
L'  in  Eg  (2.11)  =  L  in  Eg  (2.5)  (2.15) 
R'  in  Eg  (2.12)  =  [  R  in  Eg  (2.  6)]'"  (2.16) 


These  relationships  allow  the  corresponding  convergence 
types  for  H-functions  defined  by  Eg  (2.2)  to  be  written. 
Table  2  lists  these  types  of  convergence. 
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*  Table  2.  Convergence  Types  for  H-Functions  of  Eg  (2.2) 


TYPE 

D' 

E' 

L' 

H(cz) 

Izl 

|arg(z)| 

I 

>0 

<0 

>E‘if 

-1 

RHP 

res 

>0 

<  H 

■<-T^ 

II 

kO 

<0 

-z 

RHP 

res 

>0 

<  H 

III 

>0 

>0 

>E'y 

IHP 

res 

>0 

<  n 

IV 

kO 

>0 

lilP 

res 

>0 

<  n 

VI 


>0  =0  iO 


iO  =0  <0 


RHP  res 

IMP  res 

-j;  RHP  res 
IHP  res 


c 

±1 

c 

K' 

c 

R' 


<  n  ,  < 


JTD' 


<  K  ,  i 


HD' 


v^re,  if  L'  <  -1  in  Type  VI  convergence,  one  nay  use  the  sum 

■>' 

of  either  lilP  or  RHP  residues  at  |z|  = 

c 

As  expected,  there  is  a  conplete  interdiange  of  lilP  and 
RHP  poles.  Type  I  convergence  for  the  H-function  of  Eg  (2.2) 
in  Table  2  corresponds  to  Type  III  convergence  for  the  H- 
function  of  Eg  (2.1)  in  Table  1.  Similar  statements  ^pply 
between  Types  II  and  IV,  III  and  I,  and  IV  and  II.  Even  Types 
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V  and  VI  in  Table  2  correspond  to  Types  V  and  VI  in  Table  1 
fdien  considering  the  relationship  between  R  and  R'  in 
Eq  (2.16). 

Either  set  of  conditions  is  sufficient,  but  not  necessary 
for  the  appropriate  H-function  to  converge.  Still,  nearly  all 
of  the  nany  special  cases  of  the  H-function  satisfy  the 
convergence  conditions.  Furtter,  the  conditions  consistently 
and  correctly  identify  the  valid  range  of  the  variable  over 
which  the  H-fxmction  representation  equals  the  special  case. 

One  should  not  think  of  these  conditions  as  constraints  on 
the  H-fvnction,  as  Springer  [1987]  did.  Instead,  they  should 
be  viewed  as  a  sufficient  tool  to  deterndne  how  the  H-fvnction 
can  be  evaluated  with  the  svm  of  residues  -  those  in  the  LHP  or 
those  in  the  RHP.  Dependiitg  on  vdiere  Jordan's  Lemna  is 
satisfied,  residues  at  the  IHP  poles,  or  RHP  poles,  or  either 
are  svimed  for  different  values  of  the  variable.  These  choices 
are  succinctly  given  in  Table  1  (and  Table  2)  for  most  cases  of 
interest . 

To  demonstrate  that  the  convergence  conditions  are 
sufficient,  but  not  necessary,  consider  the  following 
representation  of  e^  for  x>0  as  a  Meijer  G-fvmction  [Pnidnikov 


et  al,  1990,  p.  633] 
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=  TiirW  '•  (“)'  ]  <2-i'') 

(2.18) 

where  d  is  an  arbitrary  constant.  For  this  H-f unction,  D=E=-1, 
L=  "  does  not  meet  any  of  the  types  of 

convergence  listed  in  Table  1  above.  Using  the  definition  in 
Eg  (2.1),  this  H-function  can  be  written  as 


n  1 

sin  (dji)  2ni 


r(s) 


r(i-d+s)  r(d-s) 


(x)  ®  ds 


(2.19) 


Using  the  reflection  formula  for  the  ganma  function  [Abramowitz 
and  Stegun,  1970,  p.256,  6.1.17]  with  z=d-s,  this  is 


n 


sin  (drc)  2rci 


_ 

III 

n 


r(s) 


CSC  (7l(d-s)) 


(x)“®  ds 


(2.20) 


1  _  1  si 

sin  (Ctrl)  2iii 

'  n 


sin  (nd-rrs)  r(s)  (x)  ®  ds 


(2.21) 


sin 


(sin(nd)cos(ns)  -  cos(iid)sin(iis)) 


r(s)  (X)'®  ds 


(2.22) 


Now,  r(s)  has  LHP  poles  of  order  1  at  Sj=  -J,  J=0,1, _  These 
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are  the  only  poles  of  the  integrand.  At  these  values  of  s, 
sin(ns)  vanishes  and  cos(ns)  =  Using  the  residue 

theorem  to  evaluate  the  contour  integral  as  the  sxxn  of  residues 
produces 

=  sin^(drc)  ^  [(-l)^sin(nd)j  [  ^  J;  ■] 

(2.23) 

=  e*  Q.E.D.  (2.24) 

This  proves  that  the  convergence  conditions  given  in  Table  1 
are  sufficient,  but  not  necessary,  for  the  H-f unction  to 
converge. 

2.4.  PROPHCTIES 

Among  the  useful  properties  of  the  H-function  eure  the 
identities  dealing  with  the  reciprocal  of  an  argunent,  an 
argunant  to  a  power,  and  the  multiplication  of  an  H-ftmction  by 
the  argument  to  a  power  [Carter  and  Springer,  1977]. 
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2.4.1.  RECIPROCAL  OF  AN  mSOHESTT 


[  j=i,...,q  i=l,...,p  J 

2.4.2.  ARGtHEMT  TO  A  POWER 


(2.25) 


(2.26) 


for  c  >  0 
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2.4.3.  MJLTIPLICATIOW  BY  TOE  ARGtMaiT  TO  A  POWER 


P  q 


'  =  {M]  ■- 

i“l/ • • • /P  j*!/ • • • 


p  q 


z  : 

ffa.+A.c,A.l\  ;  . 

U  1  1  ij/  ’ 

11  3  3  3}J 

i=l , . . . ,p 

j=l,...,q 

(2.27) 


2.4.4.  FIRST  REDUCTICM  PROPERTY 

If  a  pair  of  "A"  terms  axid  a  pair  of  "B"  terms  in  an  H- 
functicn  are  identical  and  one  is  in  the  nunerator  of  the 
integrand  and  the  other  is  in  the  dencndnator/  then  it  is 
equivalent  to  an  H-function  of  lower  order.  Specifically, 
[Mathai  and  Saxena,  1978,  p.  4] 


iM)  ■  '  K'*i) 

«  -issl 


.  m  n-1 
'p-1  q-1 


z 


:  {(.i-Aj)}  ;  {(b,,B.)} 

i~2, ... ,p  j—l, ... ,q“l 


(2.28) 


provided  n>0  and  g>m.  Also 
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<a 
P  q 


i=l,...,p-l  j=l,...,q 


_  jn-l  n 

■  Vi  q-1 


i=l,...,p-l  j=2,...,q 

provided  in>0  and  p>n. 


(2.29) 


2.4.5.  SEOOMD  RmUCTlOM  PROPERW 

Bodenschatz  and  Boedigheimer  [1983,  pp.  11-12]  discovered 
another  way  in  idiich  the  H-£uncticn  can  reduce  to  one  of  lower 
order,  at  least  in  the  limit.  If  any  or  B^  is  close  to 
zero,  that  ganma  term  in  the  integrand  of  Eq  (2.1)  is 
essentially  a  constant.  Thus, 


q[*l  “  rp>])  9-1 


p  q- 


i=l,...,p  j=2,...,q 


(2.30) 

for  B^  Af  0  and  m>0.  Here,  the  symbol  w  means  the  limit  of 
as  B^-^0  is  given  by  the  right  side  of  Bq  (2.30). 
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r(i-bj 


*p  q-1 


:  {[n,.,^]]  ; 

i— 1 f • • • »p 


for  B  M  0  and  nKq. 

q 


"t*)  «  r(l-aj) 


m  n-1 
p-1  q 


j=i,...,q-i 


••  {(•*'*!)}  = 
i=2,...,P 


(2.31) 


iM) 

j®lf • • • ^q 


(2.32) 


for  A.  w  0  and  n>0. 


!?>]<. - 

» ’  r  (a^ 


„  m  n| 

Vl  q 


i~if ... »p“i  • • • ^q  J 


(2.33) 


for  0  and  n<p. 


*  2.4.6.  HIERAROIICAL  WgrATTOWfiHTPS  AMOWG  ORDTO  CF 

H-PUMCTIOMS 

One  significant,  but  stxrprisingly  sinple,  discovery  was 
that  fdiole  classes  of  H-functions  are  cntedded  in  other 
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classes.  For  exanple,  the  ^  class  of  H-f unctions  is  a 
proper  svibset  of  both  the  ^  and  2  classes  of  H-f  unctions. 
These  classes,  in  turn,  are  enisedded  in  other,  hi^^r  orders  of 
H-function  classes. 

This  newly  discovered  hierarchical  relationship  among 
classes  of  H-functions  results  in  several  new  ways  in  fdiich 
certain  H-functions  can  reduce  to  H-functions  of  lower  order. 
Further,  tdiile  the  first  and  second  reduction  pri^rties  listed 
above  are  readily  apparent  and  easily  understood,  the  new 
reduction  properties  are  less  transparent. 

The  new  results  are  based  on  the  Duplication  Formula, 
Triplication  Formula,  and  Gaiiss*  Multiplication  Formula  for  the 
ganrna  function  [iUoramowitz  and  Stegun,  1970,  p.  256].  The 
Duplication  Formula  is 

-  1  2w-  ^ 

r(2w)  =  (2n)  2  r(w)  rjw+-^j  (2.34) 
which  can  be  rewritten  as 


(2.35) 

(2.36) 
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r(w)  =  2^ - ■-I——-  4"”  (2.37) 

r(4«,) 

Any  ganiiB  fxmctian  present  in  the  integrand  of  Eg  (2.1)  in 
the  definition  of  the  H-f\ncticn  can  be  replaced  with  an 
equivalent  expression  as  in  Eq  (2.35),  Eq  (2.36),  or  Eq  (2.37). 
Terns  of  the  above  equations  idiich  do  not  involve  ganna 
functions  can  be  ccnioined  with  the  parameters  k  or  c  in  the 
definition  of  the  H-function.  Using  Eq  (2.35),  a  ^  H- 

function  nay  be  rewritten  as  a  2  H-function.  Using 

Eq  (2.36)  or  Eq  (2.37),  a  ^  H-function  nay  be  written  as  a 
^  H-fmction.  The  H-functions  resulting  from  Eq  (2.36)  and 
Eq  (2.37)  appear  distinct,  but  can  be  shown  to  be  equivalent. 

I  will  state  these  new  upgrade  and  reduction  results  for 
first  order  H-functions  as  theorens  and  provide  proofs.  The 
proofs  sirply  use  the  argiment  to  a  power  property  in  Eq  (2.26) 
to  change  B  to  unity,  the  definition  of  the  H-fxmction  in 
Eq  (2.1),  and  one  of  the  forms  of  the  duplication  property  in 
Eq  (2.35)  to  Eq  (2.37). 

*  Theorem  2.1.  hJ  cr  ;  ;  (b,B)  1 


-|^  hJ  4®cz  ;  [b-^,Bj  ;  (2b-l,2B)  J 


(2.38) 
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Proof: 


cz  :  ;  (b,B)  j  =  — g-  (cz)  :  ;  (b,l)  j 


ssr  j 


S 

"B" 


■  r(bfs)  (cz)  ds 


1 

■BScT 


f  .  4-  («) 

b  4 


a 

■“B“ 


ds 


using  Bq  (2.36)  with  w  =  b-— 15-+S 


/»  1 

r 

2b-l+2s 

J 

o~ 

*  1  " 
b“  m">^3 

s 

■“B" 


ds 


=  hJ  4®cz  ;  [b-^.Bj  ;  (2b-l,2B) 


Q.E.D. 


Theorem  2.2.  ^  ^  j 


(2.39) 


2/5? 


IE"  “1  1 


4®  cz  ; 


^bf  2  J 


(2b.2B) 


(2.40) 
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Proof: 


a 


r(b+s)  (cz)  ds 


1  2. 
T2ri" 


2/i7 


r 

2(bfs)' 

r 

'  1 

b+— >^+s 

1  c 

r 

4‘®  (cz)  ds 


using  Eq  (2.37)  with  w=bfs 


4®  cz 


{2b,2B) 


(2.41) 


Q.E.D. 

Equivalence  between  Eq  (2.41)  and  Eq  (2.39)  can  be  proven 
by  xising  the  telescoping  property  of  the  gatnna  function 
r(z-fl)  =  zr(z).  The  core  of  the  proof  is  provided  below. 
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2/R 
f  T 


r(2b4-2Bs) 

bf-^+Bs 
.  « 


2^ 

4^ 


(2b-l+2Bs)  r(2b-l+2Bs) 
b--i-+B8  r  b^-^+Bs 

y  V  4 


iS 


FT 


r(2b-l+2Bs) 

- 1 - T 

r  b— ^+Bs 

.  ^  - 


(2.42) 


Theorem  2.3.  cz  :  ;  (b,B)  j 


,b-l 


^  “o  2[  ^  — r)  ] 


Proof: 


A 

hJ  cz  :  ;  (b,B)  J  *  “J  J[  (c*)”®”  =  '*  ] 


s 

nr 


=  r(Ms)  (cz)  ds 


(cz)  “  ds 


using  Eq  (2.35)  with  w  = 

'‘-I  . 

“  BSiT 

,b-l 


bfs 


rl4'4-)  (='«r  - 


-2  1  H? 


1 

.T" 
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=  ”o  2[  2  •  '  ("T — 7“) '  — r)  ] 

Q«E«Da 

upgrade  and  reducticn  results  for  the  more  gaieral  class 
of  H-functions  ^  can  be  proven  using  the  same  stqps  as  in 
the  above  proofs  by  working  on  any  gamm  function  present  in 
the  integrand  of  Eg  (2.1).  Since  the  proofs  are  very  similar 
to  those  already  provided,  Theorems  2.4  througdi  2.7  will  be 
stated  without  proof.  The  generalized  results  are 


*  Theorem  2.4.  ifj  “ 

P  9 


cz 


1*1, ... ,p  j*l/ ... ,g 


i*l, ...  ,p 


j=2,...,q 


(2.45) 
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_  2  ^  Jtrt-l  n 

i=l,...,P 

/S  P  «+• 

(^4.).  1^4.),  Ml 

j=2,...,q 

for  iA>0. 

*  Ttieoran  2.5.  if?  ” 

p  q 

« 

)} 

1*1/  •  •  •  /P  •  •  • 

'9  . 

*P 

4  n 

2/il  P  q+1 

**A 

4  '’a  :  .  (2»p'“l>) 

i=l,...»p-l 

/5?  „  m  n 

■aC=r-  Vl  q 


■  (M) 


j=i,...,q 


(2.46) 


(2.47) 


(2.48) 


for  p>n. 


Theorem  2 


•6.  =  {M}  =  {(Vd)} 

i— Xf • • • ,p  j“l» • • • 


^  ^  ^  “SJ“  *  (2V^'^)'  {(®i'^i)}  '■ 

i“2 / • • • /p 


« 

{[bj,Bj},  [•I'J-'*!) 

3“X, ...  /<! 


(2.49) 


1  „  m  n+lf 

"T—Vl  q  ^  “ 
2  V? 


i~2f ... ,p  j~l/ • • • »q 


(2.50) 


for  q>in. 

Similar  results  are  availad>le  by  using  the  triplication 
formula  or  Gaims'  multiplication  formula.  For  example,  using 
the  triplication  formula,  it  is  possible  to  show  the  ^  class 
of  H-functicns  is  also  a  proper  subset  of  both  the  ^  and 
^  classes  of  H-functions. 

*  2.4.7.  GEMERALI2IWG  OOWglMIT 

Bodenschats  and  Boedi^diieinier  [1983;  Bodenschatz  et  al. 


1990]  gave  the  first  H-ftmcticn  representaticxis  «diich 
recognized  the  existence  of  a  positive  generalizing  constant, 
u.  They  showed  that  for  certain  irathenatical  fimctions  and 
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statistical  distributions,  the  A^and  parameters  need  not 
equal  unity.  They  gave  some  new  generalized  H-function 
representations,  but  did  not  state  general  conditions  under 
whioh  a  generalizing  constant  was  possible.  These  general 
conditions  will  be  given  below  as  four  theorems. 


Theorem  2.8.  If  in>0,  p>n,  and  a^-b^=l,  then 


t/"  n 

P  q| 


i=l,...,p-l 


j=2,...,q  J 


p  q 


cz  :  iM)  ,  Jubj^+l,uBjj  ; 
i=l,...,p-l 


j=2,...,q 


(2.53) 


for  u>0. 
Proof; 


h"'  " 

"p  q 


i=l,...,p-l  j=2,...,q 
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52 
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Theorem  2.10.  If  nX),  q>in,  A.=B  ,  and  a,-b  =1,  then 

A  q  1  q 


56 


(cz)  ®  ds 
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i=2,...,p 


j— 1  /  •  •  •  /<I~1  _ 


(2.58) 


Q.E.D. 


Theorem  2.11.  If  n>0,  q>m,  B  =A, ,  and  b -a,=l,  then 

q  T.  q  1 

*  ■ 

•  •  •  /P  3^1/  •  •  •  /Q""l 

"  “T  qf  ®*  *  ('**l''*l)'  {(®i'\)}  ' 


i=2,...,p 


j~l» • • • »q“i 


(2.59) 


for  u>0. 


Proof: 


»  • 
*?q  •  (®1'^)'  {(®i'^i)}  '  '  (v^'^) 

i=2,...,p  j=l,...,q-l 
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(cs)"*  ds 


59 


_  1  ^  n 
“  u  "p  q 


cz 


i=2,...,p 


,  (ua^n-uHi) 


(2.60) 


Q.E.D. 

2.4.8.  DERIVATIVE 

It  is  well  known  that  the  derivative  of  an  H-£uncticn  is 
an  H-function  of  hi^^ier  order.  Let  the  r^  derivative  of 
1^  q[z]  be  denoted  by  If  all  real  are  such  that 

-b.  .. 

- —  <  1  for  j=l,...,in  then  the  r^  derivative  of  if™  ”[z]  is 

Bj  B  ** 

[Mathai  and  Saxena,  1978,  p.  7;  Cock,  1981,  p.  83] 
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H  [8]  - 


1=1/ • • • /P 

{(v®d-=j)}  ■  «’•"> 

j=i/...,q 


(2.61) 
th 


Cook  [1981/  p.  83]  gave  an  inproved  fozinula  for  the  r 
derivative  of  an  H-fxmction  if  any  real  is  such  that 

—  ^  1  for  j=l, . . .  ,in.  The  formula  is  based  on  the  value  of 


B 


1/  defined  as 


1  =  maxiimin 


iimin  [  0/ 


j=l/...,in|^ 


largest  integer  less  than 


-b. 


(2.62) 


The  r^  derivative  of  ^[s]  is 


=  (-1)'  1^1 


:  (-I-r,l),  •*!)}> 

i=l/ ... /P 


(-r/1)  ;  {[bj-rBj,BjJ/  (0,1) 
3=1/ ... /Q 


(2.63) 
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2.4.9.  LAPIACE  TRANSFORM 

It  is  well  Imown  that  the  Laplace  transfoxm  of  201  H- 


fimcticn  is  an  H-£tmcticxi  of  hi^^ier  order.  If  all  real  are 
stich  that - =-  <  1  for  j=l, . . .  ,in  then  the  Laplace  transform  of 


q[cz]  is  [^ringer,  1979,  p.  200;  Cook,  1981,  p.  35] 


(2.64) 

Cook  [1981,  p.  82]  gave  an  inproved  formula  for  the 
Laplace  transform  of  an  H-function  if  any  real  is  s\ich  that 

- ^  1  for  j=l,...,m. 
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H<«)  }  = 


r  :  (1,1), 


{(l-bj-Bj.Bj)}  :  (I.l),  {(l-ai-*l'»i)}’  <“'« 

j— 1,  • .  •  ,<i  i~l,  •  •  •  ,P 


(2.65) 


i^re  I  is  given  by  Eq  (2.62). 


2.4.10.  POIRIER  *naNSPQRM 

It  is  well  known  that  the  Fourier  transform  of  an  H- 
function  is  an  H*-f\]nction  of  higher  order.  If  all  real  are 

such  that - iL  <  i  for  j=l, . . .  ,m  then  the  Fourier  transfom  of 

qtca]  is  [Springer,  1979,  p.  201;  Cook,  1981,  p.  35] 


ft 

5t{  H(c*)  J  r 


e^^*  H(cz)  dz 


1  j,n+l  m 
"c"  "  q  P+1 


-  -I-  *  =  ' 

j=l,...,q 


(0.1), 

i~l, ... ,p 


(2.66) 

Cook  [1981,  p.  82]  gave  an  iirproved  formula  for  the 
Fourier  transform  of  an  H-function  if  any  real  b^  is  siich  that 
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-b^ 


~  2  1  for  j=l , . . .  ,in. 


»t{  H(«) 


-  t  :  (1,1), 


{[l-bj-Bj,BjJ  ;  (1,1),  (0,1) 

j— 1,  •  •  •  ,<i  i~l,  •  •  •  ,p 


(2.67) 


vdiere  I  is  given  by  Eg  (2.62). 


2.4.11.  MPT.T.TW  TRMISroBM 

Since  the  definition  of  the  H-functicn  in  Eg  (2.1)  nay  be 
recognised  as  the  inverse  Mellin  transform  given  by  Bg  (1.10), 
the  Mellin  transform  of  an  H-functicn  is  readily  obtained. 

r* 

H(cz)  ^  =  I  s®”^  H(cz)  ds 
'^0 


(2.68) 


2.5.  SPECIAL  CRSES  -  MMHa^TICM.  FUWCTIOIfS 

With  the  parameters  given  below,  the  H-function  can 


64 


represent  all  of  the  fol lowing  nathematical  functions  as  a 
special  case  [Mathai  and  Saxena,  1978;  Code,  1981;  Bodenschatz 
and  Boedidieimer,  1983;  Bodensdiatz  et  al  1990;  Prudnikov  et 
al,  1990;  Cook  and  Barnes,  1991].  In  some  of  the 
representations,  u>0  is  a  general iziitg  constant. 

Exponential  and  Power  Functions: 

e*  =  — r-”..  —  2C  z  :  (l-d,l)  ;  (0,1),  (l-d,l)  ] 
sin  (dn)  ^ 

e"®  =  J[  z  :  ;  (0,1)  ] 
z**  e‘*  =  hJ  J[  z  ;  ;  (b,l)  ] 


=  hJ  JC  *  :  (btl,i)  ;  (b,l)  ]  |*|  <  1 

*  =  u  z  ;  (\jb+l,u)  ;  (ub,u)  ]  |z|  <  1 

*  z*"  =  :  (bfl,l)  ;  (b,l)  ]  |z|  <  M 

*  =  u  hJ  :  (ub+l,u)  ;  (ub,u)  J  1*1  ^ 

*^  =  hJ  J[  z  ;  (bfl,l)  ;  (b,l)  ]  1*1  >  1 

*  =  u  z  :  (ub+l,u)  ;  (ub,u)  ]  |z|  >  1 
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*  a*’  =  hJ  :  (bfl,l)  ;  (b,l)  ]  !*!>«: 

*  "  “  ®1  l[  "T"  *  '  (ub,u)  j  1*1  >  ^ 

(1-z)*^  =  r(bfl)  hJ  J[  a  :  (bfl,l)  ;  (0,1)  3  1*1  <  1 

*  (M-a)^  =  r(b+l)  M*’  hJ  -g-  :  (b+1,1)  ;  (0,1)  j  |*|  <  M 

(a-1)^  =  r(bH)  hJ  J[  a  :  (bfl,l)  ;  (0,1)  ]  l=|  >  1 

*  (a-l)*’  =  r(bfl)  i[  5  (b^lA)  ;  (0,1)  ]  1*1  >  e 

a^  (l-a)"^*  s  r(a+l)  hJ  J[  a  :  (b+a+1,1)  ;  (b,l)  ]  |*|  <  1 

*  a^"  (M-a)**  *  r(a+l)  ?[  "R-  •  >  (b.D  ] 

1*1  <M 

(l+z)“^  =  hJ  h  a  :  (l-b,l)  ;  (0,1)  ] 

r(b)  - 

a^  (1+a)'*  =  ]  ■  hJ-  h  a  ;  (b-a+1,1)  ;  (b,l)  ] 

r(a)  1  1 

=  n  a  :  (0,1),  (-T'^)  ] 

a  :  (l-b,l),  ' 

(0,1),  ] 
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^  2^  *  J  (bfl,l),  (a+1,1)  ;  (a,l),  (b,l)  ] 

1*1  <  1 

2[  z  :  (bH,l),  (a+1,1)  ;  (a,l),  (b,l)  ] 

1*1  >  1 


1  _a 


“B=r  ‘ 
■Fir  * 


1*1  <  1 

1*1  >1 


=  -  z  :  (b+1,1),  (a+1,1)  ; 

(a,l),  (b,l)  ] 


a+z 


n 


1  „1  1 

■arT  "l  1 


■  K)  ■  (•.4-) 


n  >  0,  |arg  a|  <  n 


1-z' 


n 


1-z' 


nn 


(4-4r)  '  ("-ir)-  (44r)  ] 


n  >  0,  m  >  1 


.2  .  _2 


-  =  2  H?  ?[  z  :  (0,1),  (0,1),  (0,1)  ; 

1  +  z  3  3 

(0,1),  (0,1),  (0,1)  3 
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z  +  2az  cos  0  +  a 


r(4-]  i-4-]  -i- 


(0.1),  (--H-'-r-]  -  (-4-'4-)  ] 

a  >  0,  |e|  <  n 


Unit  Step  Function  and  its  Conplement: 

(1  X  k  >  0 

0  0  <  X  <  k 

*  =  “l  i[  ] 

*  ”  1[  ”Tr~  *  '  (o»u)  j 

(0  X  i  k  >  0 

1  0  <  X  <  k 

*  =  ”i  l[  “T  •  ] 


(l,u)  ;  (0,u) 


Error  Function  and  its  Conplement: 
erf  X  =  -4-  r  e-®"  do 

/R  Jq 

=  hJ  2[  ] 


*  =  hJ  2[  X  :  ] 

erfc  X  =  r  d9  =  1  -  erf  x 
/S  Jx 

*  =  2[  *  •  ] 

*  =  2[  *  •  '  H"'~H '  ] 

Inccnplete  Ganim  Function  and  its  Conplonent: 

Y(a,x)  =  f*  e"®  e““^  cl0 
''O 

=  hJ  x  :  (1,1)  ;  (a,l),  (0,1)  ] 

*  =  uhJ  jC  X  ;  (l,u)  ;  (a,l),  (0,u)  ] 

r(a,x)  =  f  e®'^  d0  =  r(a)  -  7(a,x) 

'’x 

*  =  hJ  °[  X  :  (1,1)  ;  (a,l),  (0,1)  ] 

*  =  u  hJ  '  («/!)»  (0»«)  ] 
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Trigonometric  and  Hyperbolic  Fiaictions  and  their  Inverses: 
sin  a  =  “o  2[  “T  *  '  '  (®'“H  ] 

sinh  a  =  -  hJ  2[  J  ("T-'-r) '  (®'“r)  ] 

”o  2[  “T  •  '  ] 

"o  2[  'T-  ’  '  (°'“r) '  (t — r)  ] 
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arcsin  z  = 


* 


arcsinh  z 


* 


arctan  z 


-  2[  “  =  (^'4-) '  (^'4-]  = 

{4-'4-)  ■  (®'4-)  ] 

-  >4  af  “  •  {^'4)  = 

(4-4-)  ■  (“'4-)  ] 
:-^.^2[.:<i.u),  (1.4.)  •• 

[4,4.),  (O.u)] 

4-  "2  2[  ‘  =  (^'4-) '  (4-'4-)  '■ 
(4--4-) '  (®'4-)  ] 

2^  z  :  (l,u),  ; 

(-r'4-)'  ] 


arctanh  %  - 


-  -T  2[  =  (^■4-] '  (-4'4-) ' 

(4— f) '  {®'4-)  ] 

=  -  ^2  2[  ' 

(4-4-) '  ] 

aroeot  «  =  4  ^  al  *  ’  (4'4) '  ' 

(44) '  (®4)  ] 

'442[* '  (44)'  '■ 

(4'4)'  ] 

arccoth  s  =  **2  2[  *  (4'4') '  (^'”2")  ' 

(4’4)  ■  (®'4)  ] 

'  =  4- ^  K  i*  =  (4-4)  •  ' 

e’'*  sin  (te)  =  nl^  2^  h/2  8  :  [o^-?-)  (°'"H  ] 
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sin  (az)  =  hJ  2^  ^  2  j  | »  (0»1)  | 

Logarithidc  Functicns: 

-  ®[  z  ;  (1,1),  (1,1)  ;  (0,1),  (0,1)  ] 

lnz=<  0<z^l 

hJ  z  :  (1,1),  (1,1)  ;  (0,1),  (0,1)  ] 


z  >  1 
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-u^  2^  *  •  '■  (0,u),  (0,u)  ] 

0  <  z  ^  1 

2^  *  •  (!'“)»  (^»“)  »  (0/U)»  (0,u)  ] 

z  >  1 

In  (1+z)  =  z  :  (1,1),  (1,1)  ;  (1,1),  (0,1)  ] 

=  uH^  z  :  (l,u),  (1,1)  ;  (1,1),  (0,u)  ] 

In  [l+-Y-)  =  2^  *  5  (0,1),  (1,1)  ;  (0,1),  (0,1)  ] 

=  ul^  z  :  (0,1),  (l,u)  ;  (0,u),  (0,1)  ] 

In  |l-z|  =  n  ^  z  :  (1,1),  (1,1),  ' 

(1,1),  (0,1),  J-^,lj  ] 

=  «i  3^  z  :  (l,u),  (1,1),  [-y-^l)  ? 

(1,1),  (0,u),  (■^'^)  ] 

In  jl^l  =  **  3[  *  •  (®'l)'  (I'l)'  ("^'^)  ' 

(0,1),  (0,1),  (4-)  ] 
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=  3[  *  *  ' 

(O.u),  (0.1),  i-^-]  ] 

In  Jl  +  2*  cos  e  +  =  2b  3J  z  :  (1,1),  (1,1), 

=  2»i4|[  z  :  (l,u).  (1,1), 

In  Jl+z^  =  “r  2[  *  •  (^'^)'  ' 

[1,4.) ,  (0,4.)  ] 

~  u  2^  *  *  ^  (®»'*)  j 


Bessel  Functions: 
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j“(*)  =  hJ  jC  s  ;  ;  (0,1),  (-v,u)  3 

(Maitland's  generalized  Bessel  function) 


Hypergeometric  Fmctions: 

M(a,b,-z)  =  (a;b;-z) 

=  hJ  i[  z  :  (l-a,l)  ;  (0,1),  (l-b,l)  ] 

r(a)  A  2 

(Cknfluent  Hypergeometric  fiaiction) 

2^1  '  -r^  ^  '• 

(0,1),  (1-C,1)  ] 
(Hypergeometric  fwction) 


P^ 


mn  - 


_  J/h) 


(0.»,  ] 


for  p  <  q  or  for  p=q+l  and  |z|  <  1 
(Generalized  Hypergeometric  functions) 
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(0,1),  {(l-bj-Bj}  ] 

(Maitland's  or  Wrist's  Generalized  Hypergeometric  functicn) 


MacRobert's  E-Pwcticn: 


e[p;  {aj;  q;  z]  =  J[  z  :  (1,1),  {(bj-Bj)}  : 

{(^•^)}] 


Meijer's  G-Function: 


*  2.6.  FUWCTIOWS  REPRESEWTH)  OVER  A  REgntlgraP  RSSS 

For  some  of  the  special  cases  of  the  H-ftnction  listed 
above,  the  H-fmction  represents  the  special  case  only  for 
certain  values  of  the  variable  z.  For  other  values  of  the 
variable,  the  H-function  takes  the  value  zero.  These  cases  can 
be  identified  by  a  restriction  on  the  variable  such  as  |z|  <  1. 

These  restrictions  arise  front  the  convergence  conditions 
for  the  H-function  given  earlier.  These  H-fxmctions  are  of 
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Convergence  Type  VI  vdiich  means  they  can  be  evaluated  by  the 
sum  of  IHP  residues  for  |z|  <  ^  and  by  the  negative  sum  of 
RHP  residues  for  |z|  >  ^  ♦  But  the  ^  and  2  classes  of 
H-fmcticns  have  no  RHP  poles  so  the  value  of  the  negative  sixn 
of  RHP  residues  is  zero.  Therefore,  H(z)=0  for  |z|  >  ^  . 
Similarly,  the  H^  ^  and  H^  2  classes  of  H-functions  have  no  IHP 
poles  so  the  value  of  the  sum  of  IHP  residues  is  zero. 
Therefore,  H(z)=0  for  |z|  <  . 

Through  scaling  the  variable  with  the  paraometer  c,  it  is 
possible  to  change  the  value  there  the  H-functicn  changes  from 
representing  the  special  case  to  taking  the  value  zero.  For 
example,  an  ^  H-function  can  exacstly  represent  the  power 
function  z^  for  |z|  <  M  there  M  is  any  finite  positive 
constant.  Provided  M  is  finite,  M  may  be  as  large  as  desired. 
Similarly,  an  H-function  can  exactly  r^resent  the  same 
power  ftnchion  z^  for  |z|  >  c  there  e>0  may  be  as  snail  as 
desired.  These  scaled  H-function  representations  (given 
earlier  in  the  list  of  special  cases)  allow  a  nearly  complete 
r^resentation  of  the  special  cases. 

Another  way  to  avoid  this  limitation  of  Convergence  Type 
VI  H-functions  is  to  allow  a  sli^tly  different  ftnction  to  be 
represented  for  certain  values  of  the  variable.  The  method 
basically  involves  introducing  poles  into  the  other  half  plane 
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so  that  trtien  the  H-functicn  is  evaluated  by  sunmong  the 
residues  at  these  poles,  a  function  very  close  to  the  desired 
special  case  is  detained. 

An  inportant  consideration  in  this  a^roacdi  is  to  ensure 
the  LHP  and  RHP  poles  nay  be  properly  divided  by  a  contour  C. 
I  will  give  several  examples  of  this  teohnique  for  the  power 
function  x^.  In  some  of  these  representations,  u>0  and  v>0  are 
generalizing  constants.  In  all  of  these  representations,  e>0 
may  be  as  small  as  desired.  The  first  two  examples  involve  a 
horizontal  shift  of  either  the  RHP  poles  of  ^[x]  or  the  LHP 
poles  of  nji  ^[x]  by  e. 


* 


0  <  X  <  1 
X  >  1 


c  X  ;  (bfl-e,l),  (btl,l)  ; 

(b,l),  (b-£,l)  ] 


* 


(ub+l-\K,u),  (vb+l,v)  ; 
(^,v),  (ub-uE:,u)  ] 


r 


_b+e 


* 


0  <  X  <  1 
X  >  1 


€  X  ;  (b+1,1),  (bfl+£,l)  ; 


(bfc,l),  (b,l)  ] 
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*  =  uvc  2^  X  :  (iib+l,u),  (irfb+ve+l,v)  ; 

(irt>+v£,v),  (tib,u)  ] 

The  next  two  exainples  involve  a  vertical  shift  of  either  the 
I3IP  poles  of  j^[x]  or  the  LHP  poles  of  ^[x]  by  c.  This 
introduces  both  conplex  poles  and  a  ccnplex  function  for 
certain  values  of  ^le  variable. 

•=  ie  X  :  (b+l-i£,l),  (bfl,l)  ; 

« 

(b,l),  (b-ie,l)  ] 

«  uvic  2^  ^  *  (ub4-l-uie  ,u) , 

(vb+l,v)  ;  {yb,v),  (ub-uie,u)  ] 


0  <  X  <  1 

=  ie  2^  *  '  »  (b^l+iCrl)  ; 

X  >  1 

4 

(b+ic.l),  (b,l)  ] 

=  uvie  2C  *  •  («b+l»u), 
(^d>fvie+l,v)  ;  (vb+vic,v),  (ub,u)  ] 


Although  these  restrictions  on  the  variable  z  have  always 
been  present  for  Type  VI  convergent  H-functions,  this 
limitation  of  the  H-fimction's  ability  to  represent  certain 
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fxnctions  over  all  values  of  the  variable  was  only  recently 
discovered.  The  practical  methods  siiggested  above  to  ndnindze 
the  iirpact  of  this  limitation  through  scaling,  a  slightly 
different  function,  or  use  of  conplex  parameters  are  also  newly 
developed. 


*  2.7.  MOMEMTS  OF  Ml  H-FUNCTIOM  MID  INFINITE  SatjABILlTY 

For  a  function  of  a  real  variable,  f(x),  idiich  is  nonzero 
only  for  positive  values  of  the  variable,  the  r^  moment  about 
the  origin  of  f(x)  is  given  by 

u  =  f  f(x)  dx  (2.69) 

r  Jq 

provided  the  integral  in  Eg  (2.69)  exists.  If  f(x)  can  be 
represented  as  an  H-function,  it  is  often  easier  to  find  the 


r^  mcment  of  f(x)  using  the  Mellin  transform  formula 


Mr  = 

In  particular,  if 

f(x)  =  “I  CX 


i=l,...,p  j=l,...,qj 


(2.70) 


(2.71) 


then  using  Eg  (2.68) 
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r 


1 

m 

n  r 

j=i 

n 

n  r 

i=l 

p 

n  r 

a.+A  +A.r 

q 

n  r 

K'W) 

i=n+l 

L 1  1  1  j 

j=mfl 

(2.72) 

Use  of  Eq  (2.72)  consistently  produces  the  correct  values 
for  the  moments  of  all  the  special  cases  of  the  H~ftnctian 
presented  earlier.  But  there  are  seme  fmetions  for  tdiich 
Eq  (2.72)  will  produce  a  value  for  the  r^  moment  even  idien  the 
integral  in  Eq  (2.69)  does  not  exist.  Consider,  for  exanple, 
the  sine  f motion,  sin  x.  The  xeroth  moment,  imy  be 
interpreted  as  the  signed  area  mder  the  sine  fmetion  between 

sero  and  infinity. 

00 

=  J  sin  X  dx 

=  lim  I  sin  x  dx  (2.73) 

'•o 

The  value  of  the  integral  in  Eq  (2.73)  oscillates  between  zero 
and  two  and  does  not  approach  a  limit  as  w.  Similarly,  all 
moments  of  sin  x  do  not  exist.  However,  using  the  R-fmetion 
representation  of  sin  x  and  Eq  (2.72)  to  compute  the  moments  of 
sin  x  does  produce  valid  values.  In  particular. 
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r 

“j- 

(-1)  r!  r  an  even  integer 

0  r  an  odd  integer 


(2.74) 


This  ai>parent  ccntradicticn  can  be  resolved  with  the 
concept  of  infinite  sumability  of  integrals.  A  thorough 
discussion  of  infinite  sunrability  is  b^ond  the  scope  of  this 
thesis,  but  a  brief  introduction  to  the  concept  throush 
infinite  sums  will  be  given.  Euler  began  the  sttidy  of  this 

®  i-1 

topic  by  considering  the  value  of  the  infinite  sum  Y,  (~‘^) 

i*! 


i-1 

The  limit  of  partial  sums  S  =  ^  (-1)  does  not  exist,  since 

i=l 

S  =1  if  n  is  odd  and  S  ^0  if  n  is  even.  Euler  believed  the 
n  n 


infinite  sum  should  have  the  valije 


the  average  of  the  partial  sums  o  -  -4-  7  S. 

n  “  i=i  1 


since  it  is  the  limit  of 
n 


whether  n  is  odd 


or  even. 

Since  Euler's  time,  OesAro  and  Holder  have  developed  well- 
accepted  schemes  to  find  the  values  of  infinite  suns  idtose 
partial  suns  do  not  converge  to  a  limit.  Holder's  schannas  are 
based  on  the  averages  of  the  partial  suns  (H-1),  the  averages 
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of  those  averages  (H-2),  etc.  Gesso's  schanes  (C-1,  C-2, 
etc.)  are  less  intuitive  and  sli^tly  more  ccnplicated.  The 
abroach  described  in  the  previous  paragra^  defines  both  the 
C-1  and  H-1  schemes,  iduch  are  identical  at  the  first  level. 

It  is  inportant  to  recognize  that  if  an  infinite  sun  does 
converge,  any  of  the  sunmbility  sdiemes  will  produce  the  sane 
result.  Similarly,  an  infinite  sun  vhich  is  C-m  (H-m)  sunnnable 
will  produce  the  sane  result  using  the  C-n  (H-n)  schene  idiere 
n>m. 

There  is  a  direct  analog  to  this  approach  for  integrals 
over  an  infinite  range.  An  integral  whidh  does  not  exist  nay 
still  be  suimable  under  a  sunmbility  scheme.  The  zeroth 
moment  of  the  sine  function  is  an  example  of  such  an  integral, 
under  the  C-1  or  H-1  sumability  schemes  for  integrals,  the 
integral  for  the  zeroth  moment  of  the  sine  function  takes  the 
value  unity,  the  sane  as  that  produced  with  Eg  (2.72). 

There  are  several  functions  for  fhich  Eg  (2.72)  will 
produce  a  value  for  the  r^  moment  even  vhen  the  integral  in 
Eg  (2.69)  does  not  exist.  In  this  case,  however.  Eg  (2.72) 
always  produces  the  "correct”  value  under  an  appropriate 
suimability  schene.  It  is  as  though  the  H-function  knows  about 
infinite  sunmbility  and  uses  it  correctly  when  it  is 
apprcpriate  to  do  so.  The  first  few  moments  of  some 
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trigcnometric  and  hyperbolic  functions  are  given  in  Table  3 
below.  Althousdi  the  integrals  for  these  monents  as  in 
Eq  (2.69)  do  not  exist,  these  values  are  widely  accepted  as 
"correct." 


Table  3.  Moments  of  Trigonometric  and  Hyperbolic  Functions 


Moment 

Zeroth 

First 

Second 

Third 

Fourth 

Fifth 


1 

0 

-2 

0 

24 

0 


£OS  2 
0 
-1 
0 
6 
0 

-120 


sinh  X 
-1 
0 
-2 
0 
-24 
0 


cosh  X 
0 
1 
0 
6 
0 

120 


2.8.  EVMilATIQM  OF  TOE  H-PUWCTIOM 


Like  most  contour  integrals  in  the  conplex  plane,  the  H- 
function  is  usually  evaluated  by  sunning  the  residues  at  the 
poles  of  the  integrand.  The  contour  C,  fdiich  is  sometimes 
referred  to  as  the  Bromwich  path,  is  connected  to  a  semi¬ 
circular  arc  to  create  a  Bromwich  contour,  a  closed  curve  in 
the  ccnplex  plane.  By  the  residue  theorem,  the  value  of  the 
integral  around  the  closed  Bronwidi  contour  in  the  positive 
(counter-clochwise)  direction  is  the  sum  of  the  residues  at  the 
poles  enclosed  by  the  contour.  Under  very  general  conditions. 
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the  ccntributicn  of  the  sani-circular  arc  vanishes  as  the 
radius  increases  without  bound.  In  this  case,  the  desired 
integral  over  the  Bromich  path  C  eqpjals  the  integral  aroisid 
the  closed  Brorodch  contour,  the  sun  of  residues  at  the  poles 
interior  to  the  closed  contour. 

Since  the  Bromfic^  path  C  (hjf  definition)  divides  the  UIP 
and  FHP  poles  of  the  integrand  of  the  H-function,  connecting  a 
sand-circular  arc  to  the  left  will  enclose  all  the  UIP  poles  of 
the  integrand  as  the  radius  increases  without  bound . 
Travelling  around  this  Brorodch  contour  in  the  positive 
(counter-clockwise)  direction,  we  cover  the  Bronwich  path  C  in 
the  desired  direction  from  u-ioo  to  u+ioo.  Under  the  general 
conditions  referred  to  earlier,  the  desired  integral  along  the 
Bronidch  path  C  equals  the  sun  of  residues  at  the  UIP  poles. 

Conversely,  connecting  to  C  a  send-circular  arc  to  the 
right  will  enclose  all  RHP  poles  of  the  integrand  as  the  radius 
increases  without  bound.  Travelling  arouxi  this  BromMich 
contour  in  the  positive  (counter-clockwise)  direction,  we  cover 
the  Bromwi^  path  C  in  the  t^posite  of  the  desired  direction 
from  <i>-ioD  to  u+ioo.  Changing  the  direction  travelled  around  the 
contour  from  counter-clockwise  to  clockwise  sinply  reverses  the 
sign  of  the  integral.  Under  the  general  conditions  refeixed  to 
earlier,  the  desired  integral  along  the  Brorodch  path  C  from 
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u-ioo  to  (i)+iao  equals  the  negative  sun  of  residues  at  the  RHP 
poles. 

The  convergence  conditions  given  in  Section  2.3  take  care 
of  these  details  and  indicate  idiere  the  H-function  nay  be 
evaluated  as  the  sun  of  lilP  resi<hies  or  the  negative  sun  of  RHP 
residues.  Depending  on  the  pattern  of  the  poles  of  the 
integrand,  the  process  will  result  in  a  finite  or  infinite 
series.  For  all  of  the  special  cases  listed  above,  this 
evaluation  method  will  produce  a  series  tdiich  equals  the  Taylor 
series  expansion  of  the  special  case. 

While  it  is  possible  to  verify  the  H-fimction 
representations  of  the  special  cases  in  this  manner,  it  is 
often  difficult  and  tedious  to  produce  the  series,  especially 
with  poles  of  multiple  orders.  Eldred  [1979]  wrote  a  computer 
program  to  evaluate  the  H-fxxiction  by  sunning  the  residues  at 
the  appropriate  poles.  Cook  [1981;  Cook  and  Barnes,  1981] 
improved  the  program  and  added  extra  capabilities  dealing  with 
algebraic  combinations  of  independent  H-function  variates. 
Still,  Cook's  program)  will  evaluate  a  general  H-function  over  a 
specified  range  (and  with  a  specified  interval)  of  the  real 
variable  x  idien  the  parameters  are  input. 


CSAFTSi  3 


IHE  H-FUNCTION  DISmiBUTICN 


3.1.  D^IWITIGM 

Because  the  H-£unction  in  Bq  (2.1)  can  exactly  r^resent 
the  kernel  of  many  caiman  probability  density  functions,  it  was 
natural  to  define  an  H-fimcticn  distribution  as  the  product  of 
an  H-fuicticn  and  a  ccnstant,  k,  idiich  normalizes  the  area 
under  the  H-functicn  (over  the  sqppropriate  range)  to  unity.  An 
H-functicn  variate  has  the  fol  lowing  probability  density 
function  (p.d.f.)  [Carter,  1972;  Carter  and  Springer,  1977, 
pp.  545-546;  ^ringer,  1979,  p.  200]: 


f(x)  = 


P  9 


cx 


i— 1, ... ,P  3“ if • • • f? 


X  >  0 

otherwise 

(3.1) 


In  this  case,  the  random  variable  X  is  called  an  H-function 
variate  idiich  follows  an  H-function  probability  law  or  H- 
f unction  distribution. 

Use  of  an  H-function  representation  as  in  Eq  (3.1)  for  a 
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probability  distribution  has  several  advantages.  First,  it 
tsiifies  nearly  all  cantmon  continuous  probability  distributions 
of  positive  random  variables  \jnder  one  very  general  class  of 
functions.  The  many  named  distributional  forms  idiich  arise 
naturally  in  comnon  problems  in  probability  may  be  managed  viith 
only  one  function,  the  H-function.  The  H-function  can  also 
r^resent  an  infinite  nurrber  of  other,  imnamed  distributional 
forms. 

The  H-function  also  eliminates  the  need  to  specify  the 
range  of  the  random  variable  for  idiich  the  density  is  nonzero. 
The  H-function  exactly  represents  the  desired  density  over  the 
apprc^riate  range  and  is  zero  elseidiere. 

Many  characteristics  of  a  probability  distribution  such  as 
the  monents  about  the  origin,  the  cunulative  distribution 
function,  and  the  Laplace,  Fourier,  or  Mellin  transfonn  are 
easily  found  from  the  H-function  representation.  Finally,  if 
new  random  variables  are  defined  by  algdbraic  conbinations  of 
independent  H-function  variates,  the  densities  of  the  new 
random  variables  are  easily  obtained. 

3.2.  MQMEMTS  OF  AN  H-PUNCTIOW  DISTOIEOTIOW 

For  the  real  random  variable  X  with  p.d.f.  f(x)  thich  is 
nonzero  only  for  positive  values  of  the  variable,  the  r^ 
moment  about  the  origin  of  f(x)  is  given  by 


provided  the  integral  in  Eg  (3.2)  exists.  If  the  density  f(x) 
can  be  represented  as  an  H-function,  it  is  often  easier  to  find 
the  r^^  manient  of  f(x)  using  the  Mellin  transform  fomula 

In  particular,  if 


as  in  Eg  (3.1)  then  using  Eg  (2.68)  [Carter,  1972;  Carter  and 
Springer,  1977,  pp.  546-547;  Springer,  1979,  pp.  201-202], 


m 

r  ^ 

^  f  \ 

n  r 

b.+B.+B.r 

n  r  l-a.-A.-A.r 

j=i 

L  3  3  3  J 

i=l  1  1  1  iJ 

^  ‘  p  r  1  **  r 

n  r  a.+A.+A.r  0  T  1-b  .-B  .-B 

i=n+l  I  ^  ^  ^  ^  j=infl  I  J  J 

(3.5) 

Eg  (3.5)  is  very  useful  in  finding  the  nonents  of  all  H- 
functicn  variates,  including  all  of  the  special  cases  listed  in 
Section  3.5. 

3.3.  EVAmATIOM  OF  THE  H-POWCTICH  DlglRIBDTICM  COMSTANT 


Code  [1981,  p.  109;  Cook  and  Barnes,  1991]  gave  a 
practical  method  to  detennine  the  normalizing  constant,  k,  for 
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the  H-functicn  distributicn.  If  the  integrand  of  q[c^] 
no  pole  or  zero  at  s=l,  then  since  the  zeroth  moment  of  a  valid 
statistical  distribution  must  be  unity, 

1 


k  = 


=  c 


n  rfa  n  rfi-b  -b  ] 

i=n+l  ^  ^  n=mfl  I  J 


q 
n 

j=nH-l 


using  Bq  (3.!;)  with  r^O. 


(3.6) 


The  condition  that  the  integrand  of  Ii|^  q[cx]  has  no  pole 

at  s=l  cones  from  the  need  to  evaluate  the  Mellin  transform  at 

s=l.  A  commonly  met  condition  of  nearly  all  H-function 

-b. 

distributions  vdiich  guarantees  no  pole  at  s=l  is  that  — ^  <  1 


for  j=l, . . .  ,m  and 


^-«i 


>  1  for  i=l,...,n.  For  the  first 


order  H-function  ^[cx],  Jacobs  [1986,  p.  46]  noted  that  this 
restriction  corresponds  to  the  region  in  his  (B,b)  plot  [p.  53] 
above  the  line  b=-B.  He  also  shows  [pp.  71-72]  that  first 
order  H-function  distributions  in  this  region  are  uniquely 
determined  by  their  moments. 
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The  conditian  that  the  integrand  of  zero 

at  s=l  arises  from  the  need  to  have  a  nonzero  Mellin  transfonn 
at  s=l.  If  the  Mellin  transfonn  were  zero  at  s=l,  no  constant 
would  exist  to  create  a  valid  p.d.f . 

3.4.  OMUIATIVE  DISTRIBUTION  roNCTIOtl 

Ihe  ciznilative  distribution  function  (c.d.f.)  of  an  H- 

function  distribution  was  available  as  sinply  one  (unity)  minus 

another  H-function  of  higher  order  [Eldred,  1979,  pp.  139-140; 

Springer,  1979,  p.  243;  Cook,  1981,  p.  103].  If 

f(x)  =  Id^  ^[cx]  is  the  p.d.f.  of  the  H-function  variate  X  and 
P  <I 

F(x)  represents  the  c.d.f.  then 

F(x) '  1  -  4  ' 

(3.7) 

Code  [1981,  p.  103]  gave  another  equivalent  representation 
tdiich  allows  simultaneous  coninitation  of  the  p.d.f.  and  the 
c.d.f.  of  an  H-f\sicticn  distribution  by  the  sun  of  residues. 
'Hie  H-fxniction  in  this  representation  of  the  c.d.f.  has  a 
nearly  identical  pattern  of  poles  and  residues  as  the  H- 
function  representing  the  p.d.f. 
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F{x)  =  1  -  k  X  cx  :  (0,1)  ; 

(-l,l),{(b.,B.]}] 

(3.8) 

Cook  [1981,  p.  104]  also  showed  that  the  c.d.f.  (in 
addition  to  the  corrplementary  c.d.f.)  was  an  H-function.  In 
particular, 


P(x)  = 


•?+!  qtl[  “  •  '■ 

Ci  qn[  “  =  '■ 

(0'i)'{(V®3'Bj)}]- 


1=0 


-k 

c 


I>0 


(3.9) 

where  I  is  as  given  in  Eq  (2.62). 

It  is  possible  to  introduce  an  arbitrary,  positive 
generalizing  constant,  u,  to  Eq  (3.9).  The  new  generalized 
formula  for  the  H-function  representation  of  the  c.d.f.  of  an 


H-function  variate  is 
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e 


(3.10) 


t^re  I  is  as  given  in  Eq  (2.62)  and  u>0  is  arbitrary. 


3.5.  SPECIAL  Ca^SES  -  STRTISTICiM.  DlgntlBDTIOMS 

With  the  parameters  given  belcw,  the  H**function 
distributioD  can  r^resent  all  of  the  following  probability 
density  functions  as  a  special  case  [Carter,  1972;  Carter  and 
Springer,  1977,  pp.  547-549;  Mathai  and  Saxena,  1978,  pp.  10- 
12;  Eldred,  1979,  pp.  104-108;  ^ringer,  1979,  pp.  202-207; 
Cook,  1981,  pp.  85-87;  Cook  and  Barnes,  1981,  p.  300; 
Bodenschatz  and  Boedigheimer,  1983,  pp.  17-22;  Bodenschatz  et 
al,  1990].  The  H-f unction  representations  of  the  cinulative 
distribution  functions  and  a  formula  for  the  moments  about  the 
origin  are  also  presented.  In  some  of  the  representations,  v 
denotes  an  arbitrary  positive  constant.  In  the  c.d.f. 
representations,  u  denotes  an  arbitrary  positive  constant. 
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Gama  distributicn 
.r 


f(x|r,X}  = 


r(r) 

-A—  hJ  ?[  XX  :  ;  (r-1,1)  ] 
r(r)  ” 


X  >  0 
r,X  >  0 


F(x|rA)  =  hJ  XX  :  (1,1)  ;  (r,l),  (0,1)  ] 

=  hJ  Xx  :  (l,u)  ;  (r,l),  (0,u)  ] 

r(r)  2 

1  r(r+r*) 


X  >  0 


p  = - , - 

r(r) 


Ej^Kxiential  distribution  (Ganna  distribution  with  r=l) 

(Neibull  distribution  with  ^=l) 

£(x|X)  =  X 

=  X  hJ  J[  Xx  ;  ;  (0,1)  ]  x  >  0 

X  >  0 

*  P(x|X)  =  hJ  Xx  :  (1,1)  ;  (1,1),  (0,1)  ] 

*  =  u  hJ  jt  '  (1A)»  (0»u)  ]  *  >  0 
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=  -V  r(i+r) 

^  X 


Chi-Square  distributicn 


(Ganna  distribution  viith  r 


and  X  —  — ) 


e 


X  >  0 


u  >  0 


*  F(x|p)  =  '  j 

1 

w 

-■^al 

-5-  •  (i»i)  ; 

» 

[-3--^] 

,  (0,1)  ] 

u 

T 

L  J 

1“  ®i  2 

[  ^  :  (l,u)  ; 

1 ,  (0,u)  J 

X  >  0 
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Neifaull  distribution 

fix\p,\)  =  ^  X  e 


-  X 


F(x^^,X)  =  hJ  ^ 


(t- 

Jcr]  fi  1  11 

«0  1 

X  X  ;  ;  1^1- 

•  « 

X 

PxX 

1  1 

(■t) 

t 

.  f  1  ^ 

1  2 

X 

X  :  (1,1)  ;  |lx-p-Jx  (0,1) 

m 


•  (1»U)  J  ^1,— p— ^  ,  (OrU) 


X  >  0 


= 


r 

“T 


-f ) 


Raylei^  distribution  (fi^bull  distribution  with  ^=2) 

-Xx^ 


f(x|X)  =  2X  X  e 

=  /*  ^  J[  A  X  !  :  (-r.-!-)  ] 


X  >  0 
X  >  0 
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F(x|X)  =  hJ  /T  X  :  (1,1)  ;  (0,1)  ] 

=  u  hJ  2[  A  *  :  (i,«)  ;  (^'-r) '  ^  ° 


Kaxwell  distribution 
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(a,l),  (0,1)  ] 
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=  u  r(ygl.  gi  (a+^,1)  ; 

r(a)  ^  2 

(a,l),  (0,u)  ] 
X  >  0 

.  r(a»r)  T(a+p) 

^  r(a)  r(a+^+r) 

Power  Fmc±iGn  distribution  (Beta  distribution  tilth  ^=l) 
£(x|a)  =  a  0  <  x  <  1 

=  o  hJ  ®[  X  :  (a,l)  ;  (o-l,l)  ] 

=  vtt  hJ  J[  X  :  (v(o-l)+l,v)  ;  (v(a-l),v)  ] 

o  >  0 

F(x|a)  s  a  X  :  (1,1),  (a+1,1)  ;  (a,l),  (0,1)  ] 

-  uwx  2^  X  :  (l,u),  (UB+I,v)  ;  (w«,v),  (0,u)  ] 

X  >  0 

=  -sV 


Onifonn  distribution  (Beta  distribution  with  a=^=l) 

(Power  Function  distribution  with  a=l) 


f(x)  =  1 


=  hJ  ®[  X  ;  (1,1)  ;  (0,1)  ] 


0  <  X  <  1 
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=  V  hJ  J[  X  :  (l,v)  ;  (0,v)  ] 

F(x)  =  X  ;  (1,1),  (2,1)  ;  (1,1),  (0,1)  ] 

=  uv  jC  X  :  (l,u),  (v+l,v)  ;  (v,v),  (0,u)  ] 


X  >  0 

1 

■  T+r 

Pareto  distribution 

f(x|a)  =  a  X  >  1 

=  o  X  :  ;  (-«-l,l)  ) 

=  vtt  hJ  J[  X  :  (l-v(a+l),v)  ;  (-v(a+l),v)  ] 


o  >  0 


F(x|o)  =  a  X  :  (1,1),  (l-u,l)  ;  (-a,l),  (0,1)  ] 

=  uwx  ^  2^  X  :  (l,u),  (l-UK,v)  ;  (-w,v),  (0,u)] 


r 


a 

o-r 


for  a>r 


Half -Cauchy  distribution 


X  >  0 
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* 


* 


'  i[  -I-  =  (“--T-) '  (o-^)  ]  *  >  ® 

0  >  0 

Wx|9)  =  ^  ^[  -|-  :  (1,1),  (-J---T-)  '■ 

(4,4.),  (0,1)  ] 

=  >4 "[  -f  =  (-f '4-)  ! 


do  not  exist 


X  >  0 

if  r  is  even 

if  r  is  odd 


Half -Student  distribution 

2  r 

f(x|p)  =  - 


^  '•(4.)  (1 .  -^) 
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X  >  0 
u  >  0 


* 


X  >  0 


if  p-r^-2J  for  J=0,l,2,.... 


do  not  exist  if  p-r=-2J  for  0=0,1^ 2 
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F  distribution 


f(x|t>,u) 


)  ^ 

{4-: 

). 

(4)  . 

4-] 

1 

4-) 

(l^X-H)) 

[-^) 

•  '■w 


(4  -^■*)  ] 

X  >  0 
P,<|)  >  0 


F(x|p,u)  = 


- - - "2  ll  -5^  * 

Jl-  ;  (-^.1)*  (0,1)  ] 


u 


W  '(4) 


(i,«. 


(1-  ^,1)  ;  ,  (o.u)  ] 


X  >  0 
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r|  -  1  f  IlsT  ") 


rr  *  rr 


= 


if  -t^-3  for  J=0,l,2,... 


do  not  exist 


if  -2-  -r=>J  for  J=0,l,2,... 


Beta  distribution  of  the  second  kind 

.  “  r(a+^) 

=  (-5-1  - 

r(o)  ro) 


a  r(a)  ro) 


i[  ^  ] 


X  >  0 


a,^  >  0 


=  rw  r(p--  2[  ^  '• 

(a,l).  (0,1)  ] 


r(a)  ro) 


(I,u),  (l-^,l)  ; 


(a,l),  (0,u)  ] 


X  >  0 
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- - — - - r(a+r)  rO-r) 

r(o)  ro) 

=  ■  if  p-r^-J  for  J=0,l,2, . . . . 

do  not  exist  if  p-r=-J  for  J=0,l,2,.... 

k 

General  Hypergeonetric  distributicn 


c 


X  >  0 
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X  >  0 
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It  should  be  noted  that  since  the  H-function  distributicn 
is  oily  defined  over  positi^  values,  symnetric  and  doubly 
infinite  distributions  like  the  normal  and  Student's  t  nust  be 
manipulated  in  their  folded  forms.  While  the  normal.  Student's 
t,  and  Cauchy  distributions  are  not  special  cases  of  the  H- 
function,  their  folded  forms,  the  half-normal,  half- 
Student's  t,  and  half-Cauchy,  are  representable  as  H-function 
distributions. 

The  Venn  diagram  in  Figure  1  shows  the  relationship 
between  many  comnnon  first  and  second  order  H-function 
distributions.  The  Erlang  distribution  is  simply  a  ganma 
distribution  with  an  integer  shape  parameter  r.  The 
e^^Kxiential  distribution  is  a  proper  subset  of  both  the  Weibull 
distribution  and  the  Erlang  distribution,  vhich  is  a  proper 
subset  of  the  ganoma  distribution.  One  esqponential  distribution 
(with  X=  )  is  a  special  case  of  the  Chi-Square  distribution, 
tdiich  is  also  a  prc^r  subset  of  the  gamma  distribution.  If  a 
Chi-Square  distribution  has  an  even  nunnber  of  degrees  of 
freedom,  u,  it  is  also  an  Erlang  distribution  with  X=-^.  The 
Raylei^  distribution  is  a  proper  subset  of  the  Weibull 
distribution.  All  of  these  named  distributions,  plus  the  Half- 
Normal  and  Maxwell  distributions,  are  special  cases  of  the 
first  order  H^  ^  class  of  H-function  distributions. 
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It  is  well  known  that  the  Unifom  distribution  is  a 
special  case  of  the  Power  Function  distribution,  is  a 
special  case  of  the  Beta  distribution.  All  of  these 
distributioans  are  represented  as  ^  H-f unction  distributions. 
Further,  it  was  shown  in  Section  2.4.6  that  the  ^  class  of 
H-functions  is  a  proper  stibset  of  the  ^  class  of  H- 
functions.  Therefore,  all  of  the  named  distributions  in 
Figure  1  can  be  exactly  r^resented  as  H-function  distributions 
in  the  ^  class. 

An  H-function  r^resentation  has  not  been  given  for  the 
p.d.f.  of  the  Lognormal  or  Logistic  distributions.  Conversely, 
no  one  has  proven  that  these  distributions  cannot  be 
represented  as  H-functions.  This  is  an  area  for  future 
research. 

*3.6.  ARBITRARY  RMIGES  PC3R  TYPE  VI  H-tUNCTIOM  VARIATES 

The  H-functions  in  some  of  the  p.d.f.  representations  in 
Section  3.5  are  Type  VI  convergent  according  to  Table  1  in 
Section  2.3.  Consequently,  they  represent  the  desired  p.d.f. 
over  a  certain  range  of  the  variable  and,  since  they  lack  poles 
in  the  opposite  half-plane,  are  sero  for  other  values  of  the 
variable. 

Specifically,  the  H-ftziction  representations  of  the  Beta 
(first  kind).  Power  Function,  and  Uniform  distributions  above 
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are  all  ^  Type  VI  convergent  H-ftnctions  with  no  RHP  poles. 
As  listed  above,  they  represent  the  desired  p.d.f.  over 
0  <  X  <  1  and  are  zero  for  x  >  1. 

Similarly,  the  H-function  representation  of  the  Pareto 
distribution  above  is  an  ^  Type  VI  convergent  H-function 
with  no  lilP  poles.  As  listed,  it  represents  the  Pareto  p.d.f. 
for  X  >  1  and  is  zero  over  0  <  x  <  1. 

This  limitation  (as  discussed  in  Section  2.6)  of  the  H- 
function  to  represent  a  general  power  function  or  beta- type 
function  for  all  x  >  0  is  actually  an  advantage  Uhesi 
representing  the  statistical  distributions.  Using  the  H- 
function  r^resentation  eliminates  the  need  to  specify  the 
range  of  the  variable  for  tdiich  the  density  is  nonzero.  The  H- 
finction  representation  gives  the  desired  p.d.f.  for  the 
appropriate  range  of  the  variable  and  is  automatically  zero 
otherwise. 

The  newly  discovered  technique  of  scaling  these  Type  VI 
convergent  H-functions  discussed  in  Section  2.6  is  also 
applicable  here.  By  changing  the  value  of  c  in  the  definition 
of  the  H-fimction,  it  is  possible  to  alter  the  point  at  idiich 
the  H-function  changes  from  representing  the  special  case  to 
taking  the  value  zero.  For  the  statistical  distributions 
mentioned  earlier,  this  allows  more  flexibility  in  representing 
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the  same  functional  form  over  a  more  general  range.  Given 
below  are  the  new  H-function  representations  of  the  densities 
for  the  Three-Parameter  Beta  distribution,  the  Power  Function 
distribution  over  (0,M),  the  Uniform  distribution  over  (0,M), 
and  the  Pareto  distribution  over  (c,a>).  In  some  of  the 
representations,  v  denotes  an  arbitrary  positive  constant. 
Three-Parameter  Beta  distribution  of  the  first  kind 

f(x|a,P,M)  =  —^B=r 

r(a)  ro)  if 1 

0  <  X  <  M 

*  °  '^  ?[  4  =  ••  <“-!'»  ] 

a,^,M  >  0 


Power  Function  distribution  over  (0,M} 

(Three-Parameter  Beta  distribution  with  ^=l) 


f  (x|a,M)  =  a 


0  <  X  <  H 


■  "B"  l[  “H"  *'  ] 

=  hJ  -g-  ;  (v(o-l)+l,v)  ;  (v{a-l),v)  j 


a,M  >  0 


Unifonn  distribution  over  (0,M) 


(Three-Parameter  Beta  distribution  with  a=P-l) 

(Power  Func±ion  distribution  over  (0,M)  with  a=l) 

£(x|M)  —  1^  0  <  X  <  M 

*  =  hJ  -J-  :  (1,1)  ;  (0,1)  J 

*  =  -g-  J[  -J-  J  (l»v)  ;  (0,v)  J  M  >  0 

Pareto  distribution  over  (e,oo) 

f(x|a,e)  =  a  e®  X  x  >  c 

*  =  -f  ^  i[  -r  =  ] 

*  i[-r  =  ;  (-v(on),v)l 

a,e  >  0 

3.7.  TWWSPOBMKnOWS  CF  IHEEPHtEPIT  H-TOMCTIOW  VaRIABSS 

A  significant  advantage  of  using  the  H-fixiction 
representations  of  statistical  distributions  is  that  they  make 
finding  the  distribution  of  an  algdbraic  confcination  of 
independent  random  variables  much  easier,  carter  [1972;  Carter 
and  l^ringer,  1977,  pp.  549-557;  Springer,  1979,  pp.  207-219] 
showed  that  the  product,  quotient,  or  power  of  independent  H- 
fuiction  variates  themselves  had  an  H-function  distribution. 
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In  fact,  the  resulting  density  of  the  new  randan  variable  can 
imnediately  be  written  (as  an  H-finctiori)  using  only  the 
parameters  of  the  H-functions  in  the  algdt>raic  conbinaticn. 

^ringer  [1979,  pp.  217-219]  described  a  method  for 
finding  the  H-functicn  distribution  of  an  algebraic  combination 
(including  constants,  products,  quotients,  and  powers)  of 
independent  H-function  variates.  Cook  [1981,  p.  92]  combined 
the  three  theorene  of  Carter  [1972;  Carter  and  Springer,  1977, 
pp.  549-557;  Springer,  1979,  pp.  207-217]  into  one  very 
complicated  theorem.  Using  Code's  theorem,  the  H-f\inction 

V 

distribution  of  Y  =  n  X .  is  imnediately  available,  there  X ., 

jsl  D  J 

jsl,...,V  are  mutually  independent  H-f\nction  variates  and  P^, 
j=l,...,V  may  be  positive  or  negative.  Carter's  separate 
results  are  given  below,  not  Code's  combined  results. 
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m  tL 
h  \  q. 


fi(Xi)  = 


3=1,..., j=l,...,q^ 


X,  >0 

0  othervise 

V 

(3.11) 

for  i=l,2,...,H,  then  the  p.d.f.  of  the  randan  variable 


y  =  n  X .  is  given  by 
i=l  ^ 


R  R 


n  k.  H 

i=l  R  R 


i?i"‘  i?i"‘  r  f " 


A'*' 


(ii  'i)  '■ 


fyCy)  = 


y  >  0 


otherwise 


(3.12) 


tdiere  the  sequence  of  parameters 


3=1, 2,..., 

followed  by 


for  i=l,2,...,N 


j=n^+l,n^+2, . 


for  i=l,2, . . .  ,N 
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and  the  sequence  of  parameters 


j=l#2, . . . ,n^ 
followed  by 
j=m^+l,n^+2, . . . ,q\ 


for  i=l/2, . . . ,N 


for  1=1,2, . . . ,N 


In  effect,  the  fotnula  retains  the  gamna  terns  in  their 
previous  place  (nunnerator  or  denominator)  in  the  integrand  of 
the  H-function. 


3.7.2,  DISTRIHOTIOW  OF  A  OtIOTIEMT 

If  aixi  independent  H-function  variates  with 

densities  and  f2(x2),  respectively,  where 


fi(Xi)= 


"V  "i 

1  \ 


=i*i  =  '• 

3=1, ... ,p£  3“ i» • ■ • ,q^ 


Xj  >0 

otherwise 


(3.13) 


for  i=l,2  then  the  p.d.f .  of  the  random  variable  Y  = 


IS 


given  by 
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fyCy)  = 


k  k 


m.4TU  ni-Hiu 


=  {(a.D)}  ; 

(M)  ] 

y  >  0 
otherwise 


(3.14) 


t^re  the  sequence  of  parameters  {(-)}  is 

I  (*11'  ^l)'  ‘  '(®ln^' 

[^■^21  J  , . . . ,  (l-‘>2n^  ■2B2n^  /  , 

(®lnj+l'  ^+l),...,(*lp^'  ^p^)' 

(^"  ‘*2n^+l  '“211^+1'  ®2n^+l)  '  *  '  *’2q2  ■“2q2  '  ®2  J  | 

and  the  sequence  of  parameters  ^^e,Ej^  is 
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I  (^’ll'  ®ll)'  *  '(**lin^' 

(l-aji  -2^21 »  Aji]  /  •  •  • » (^"*2112  ' 

(^Im^+l'  ®lin^+l)  , . . . ,  '  ^Iq^)  ' 

[^~  *202+1  *^^202+1'  ^+l)'  **'(^‘  *2P2  ^2)  ) 

3.7.3.  DlgHtlBOTlOM  OF  A  VARIATE  TO  A  PCWTO 

If  X  is  an  H-function  variate  with  density  f  (x)  idiere 


k  hJJ 
p  q 

•  « 
cx  ;  ;  {(‘'j'®j)} 

1*1/ • • . /P  3*1/ • • • /q 
■  • 

y  0  otherwise 

(3.15) 

P 

then  the  p.d.f .  of  the  randan  variable  Y  =  X  is  given  by 
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fy(y)  = 


P-1  m  n 

k  c  H 

P  <1 


y  :  : 

i=l,...,P 

j=l, . . . ,q  J 


y  >  0 

otherwise 

(3.16) 


vAien  P  >  0  and 


fy(y)  = 


P-1  n  tn 
q  P 


j=l,...,q 

{{l-.^*AjP-Ai.-A.p)} 

i=l,...,p  J 


y  >  0 


otherwise 

(3.17) 


fdien  P  <  0. 


3.7.4. 


OSE  OE 


IS' 


(B.b)  PLOT  JN  FINDING  PCWHtS 


OF  FIRST  ORDER  H-PPHCTIOM  VMtlATES 


Even  first  order  H-f\incticns  (possessing  only  one  ganna 
term)  can  represent  a  wide  variety  of  distributional  forms. 
For  exanple,  an  H-fimction  distribution  can  have  a  sha^  idiich 
is  not  quite  that  of  either  a  Neibull  density  or  a  gamna 
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density.  Die  plot  given  in  Figure  2  of  the  (B,b)  paraometer 
space  for  first  order  H-functicn  distributions  by  Jacobs  [1986, 
p.  53;  Jacobs  et  al,  1987,  p.  134]  shows  the  ability  of  first 
order  H-function  variates  to  model  an  infinite  niarber  of  other, 
unnamed  distributional  fonns. 


b 
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Legend 

C  -  Chi-Square 
E  -  Exponential 
G  -  Ganna 
H  -  Half-Nomal 
M  -  Maxwell 
R  -  Rayleigh 
H  -  Heibull 
Blank  space  -  an 
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Figure  2.  Classical  Statistical  Distributions  as 
First  Order  H-Functions  in  (B.b)  Space 


Jacobs  [1986,  p.  55;  Jacobs  et  al,  1987,  p.  144]  showed 
how  the  (B,b)  plane  can  be  used  to  determine  the  H-fimction 
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distribution  of  positive  powers  of  H-f unction  variates.  Since 
Eq  (3.16)  in  Section  3.7.3  for  P>0  does  not  change  the  values 
of  m,  n,  p,  or  q,  positive  powers  of  first  order  H-function 
distributions  are  also  first  order  H-function  distributions. 
Further,  Jacobs  noted  that  all  positive  powers  of  the  H- 
function  variate  X  with  parameters  (b*,B*)  lie  on  the  line  in 
the  (B,b)  plane  parallel  to  the  line  representing  the  Neibull 
family  of  prdbability  density  functions  and  through  the  point 
(B*,b*). 

Jacobs  [1986,  p.  55;  Jacdbs  et  al,  1987,  p.  144]  gave 
three  exanples  of  the  vise  of  the  (B,b)  plane  in  finding  the 
distribution  of  the  power  of  H-fvnction  variates.  He  made  a 
slight  error  Utiea.  he  stated  that  the  square  of  a  half-normal 
random  variable  has  a  Chi-Square  distribution.  Instead,  the 
square  of  a  standard  (i.e.  o=l)  half-nonral  random  variable  has 
a  Chi-Square  distribution  with  i>=l.  He  correctly  stated  that 
the  square  of  a  Rayleigh  random  variable  has  an  exponential 
distribution  and  that  any  positive  power  of  a  Neibull  random 
variable  has  another,  different  Neibull  distribution. 

*  In  particular,  if  X  has  a  Raylei^  distribution  with 
parameter  X,  then  Y=X^  has  an  exponential  distribution  with 
parameter  X.  If  X  has  a  Neibull  distribution  with  parameters  ^ 
and  X,  then  Y=X^  tdiere  P>0  has  a  Neibull  distribution  with 
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parameters  and  X.  I  do  not  believe  these  relationships 
were  conmcnly  known. 

*  There  are  other  new  relationships  between  statistical 
distributions  available  from  \ising  the  (B,b)  plane  of  first 
order  H-ftnction  distributions.  If  X  has  a  half-nornal 
distribution  with  parameter  a,  then  Y=X^  has  a  ganma 
distribution  with  parameters  r  =  and  X  =  — If  X  has  a 

Maxwell  distribution  with  parameter  6,  then  Y=X^  has  a  ganma 

3  1 

distribution  with  parameters  r  =  and  X  =  — If  X  has  a 

Maxwell  distribution  with  parameter  6-/7,  then  Y=X^  has  a 
Chi-Square  distribution  with  parameter  p=3.  If  X  has  a 
Raylei^  distribution  with  parameter  X  -  then  Y^X^  has  a 

Chi-Square  distribution  with  parameter  p=2.  I  do  not  believe 
these  relationships  are  comnnonly  known,  either. 

*  The  (B,b)  plane  is  useful  for  finding  the  parameters  for 

powers  of  first  order  H-futKd:ion  variates.  If  P>1,  the 

parameters  for  Y=X^  can  be  found  by  traveling  down  an  imaginary 
line  toward  the  lower  ri^t  of  the  graph  and  parallel  to  the 
line  representing  the  Heibull  faonnily  of  densities.  If  0<P<1, 
the  parameters  for  Y=X^  can  be  found  by  traveling  up  an 
iiragizmry  line  toward  the  upper  left  of  the  gra^  and  parallel 
to  the  line  representing  the  Neibull  family  of  densities.  In 
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either  case,  the  parameter  B  for  Y  is  sinply  the  product  of  P 
and  the  parameter  B  for  X. 

3.7.5.  DISTRIHOTICM  OF  A  Sm 

Although  the  products,  quotients,  and  powers  of 
ind^>endent  H-f unction  variates  have  H-function  distributions, 
this  characteristic  does  not,  in  general,  extend  to  suns. 
Analogs  to  Carter's  results  in  Sections  3.7.1,  3.7.2,  and  3.7.3 
do  not  currently  exist  for  suns  of  independent  H-function 
variates. 

Becaxise  the  Laplace  transform  of  the  p.d.f.  of  an  H- 
fuxction  variate  is  readily  obtained  as  another  H-function,  it 
was  hoped  the  distribution  of  a  general  sun  of  H-function 

variates  had  an  H-function  distribution.  The  Laplace  transform 
of  a  density  is  nonzero  only  over  positive  values  can  be 

used  like  the  moment  generating  function  or  (diaracteristic 

function.  The  Laplace  transform  of  the  p.d.f.  of  the  sun  of 

two  independent  random  variables  is  the  product  of  the  Laplace 
transforms  of  the  individual  densities.  If  this  product  were 
available  as  another  H-function,  it  could  be  inverted 

analytically,  yielding  the  desired  density  expressed  as  an  H- 
function.  The  problem  reduced  to  determining  whether  the 
product  of  two  H-functions  was,  in  general,  another  H-function. 

Because  the  H-function  can  exactly  represent  nearly  every 
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cannon  mathenatical  function  and  statistical  density,  there  was 
anple  reason  to  suspect  that  the  product  of  two  H-functions 
was,  in  general,  another  H-function.  Indeed,  there  are  nany 
cases  vd^re  two  individual  functions  and  their  product  are  all 
special  cases  of  the  H-function.  Also,  the  sum  of  certain 
independent  H-function  variates  does  have  an  H-function 
distribution.  Examples  of  this  are  independent  exponentially 
or  gamna  distributed  variates  with  a  comnon  value  of  X. 

Unfortunately,  the  product  of  two  H-functions  might  not, 
in  general,  be  another  H-function.  Recently,  Prudnikov  et  al 
[1990,  p.  354]  gave  the  following  result 


j® 


jO-l  ^  t 
"u  V 


«  :  {(o.  .cj}  : 

i=l,...,u  j=l,...,v 


1—1, ...  ,p  j— 1/ • •  • 


dx 
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under  certain  general  conditions.  Conparing  the  left  side  of 
Bq  (3.18)  to  Bq  (1.7),  the  Hellin  transform  of  the  product  of 
two  H-f unctions  is  another  H'function.  Using  a  result  from 
Erd^lyi  [1954,  p.  308,  No.  13],  it  was  possible  to  verify  the 
result  in  Bq  (3.18).  Still,  Bq  (3.18)  is  not  entirely 
satisfying.  The  H-function  on  the  ri^t  side  of  Bq  (3.18)  may 
not  even  be  a  valid  H-function  or  Mellin-Bames  integral.  The 
variable  is  a,  the  transform  variable,  tdiioh  does  not  appear  in 
the  argument  location,  but  in  the  parameters  of  the  H-function. 
The  variable  a  will  appear  inside  certain  products  and 
quotients  of  ganma  functions  in  the  integrand  of  the  H- 
function. 

It  has  not  been  shown  that  the  Mellin  transform  of  an  H- 
f unction  is  another  H-function.  However,  by  E<1  (2.68),  the 
Mellin  transform  of  an  H-functicn  is  basically  products  and 
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quotients  of  ganna  functions.  Considering  Eq  (3.18),  if  the 
product  of  the  two  H-functions  inside  the  integral  produced  aoi 
H-function,  then  the  Mellin  transform  integral  would  yield  its 
products  and  quotients  of  garnna  functions.  Therefore,  the  H- 
function  on  the  right  of  Eq  (3.18)  would  have  to  represent 
those  products  and  quotients  of  ganma  functions. 

An  H-function  representation  has  not  been  given  for  a 
single  ganna  function.  An  H-function  representation  of 
products  and  quotients  of  ganna  functions  seena  even  less 
likely. 

Although  this  discussion  is  not  a  rigorous  nathenatical 
proof,  it  leads  to  the  conjecture  that  the  product  of  two  H- 
functions  is  not,  in  general,  another  H-function.  Further 
research  is  needed  to  prove  or  disprove  this  conjecture. 

Because  the  product  of  two  H-functions  was  not  available 
as  another  H-function,  a  pure  analytical  solution  for  the 
density  of  the  sum  of  independent  H-function  variates  as 
another  H-function  was  not  possible.  Analogs  to  Carter's 
results  for  the  sum  of  independent  H-function  variates  do  not 
currently  exist.  Therefore,  development  of  a  practical 
technique  to  find  the  H-function  distribution  for  the  sun  (or  a 
close  approximation  to  it)  was  necessary.  A  description  of 
this  technique  is  presented  in  the  next  chapter. 


anPTBt  4 


FINDING  AN  H-FUNCTIGN  DISTRIBUTION  FOR  THE  SQM 

Because  the  product  of  two  H-functions  was  not  available 
as  another  H-function,  it  was  necessary  to  develop  a  method  to 
find  a  close  approximation  for  the  density  of  the  stn  of  two  or 
more  independent  H-function  variates.  Ideally,  this 
approximation  would  also  have  an  H-function  representation. 
This  would  allow  the  inclusion  of  the  sum  of  random  variables 
in  more  complicated  alg^raic  combiimtions  with  other 
independent  random  variables. 

For  exanple,  suE’Pose  the  density  of  the  random  variable  Z 
is  desired,  where 


^2^3 

(4.1) 

\  =  VV*3 

(4.2) 

^2  = 

(4.3) 

^3  = 

(4.4) 

and  the  are  mutually  independent  H-function  variates  for 
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i=l,...,9.  If  the  (i=l,2,3)  were  available  as  H-functicn 
variates,  Caurter's  resxilts  in  Sections  3.7.1  through  3.7.3 
could  be  used.  Eq  (3.16)  for  a  variate  to  a  power  in 
Section  3.7.3  could  be  used  to  find  the  H-f unction  distribution 

of  .  Similarly,  Bq  (3.12)  for  a  product  of  H-function 
variates  in  Section  3.7.1  could  be  used  to  find  the  H-f  unction 
distribution  of  ^2  Y^*  Finally,  Bq  (3.14)  for  a  quotient  of  H- 
fmction  variates  in  Section  3.7.2  could  be  used  to  find  the  H- 
f unction  distribution  for  the  density  of  Z.  Even  a  conplicated 
algdnraic  combination  of  independent  random  variables  such  as 
this  becomes  almost  trivial  if  the  Y^  are  available  as  H- 
f unction  variates. 

Cook  [1981;  Cook  and  Barnes,  1981]  developed  a  method  for 
finding  the  p.d.f.  and  c.d.f.  (in  tabular  form)  of  an  algebraic 
conciliation  involving  products,  quotients,  powers,  and  sure  of 
independent  H-function  variates.  Cook  also  developed  a  FCRTBUH 
computer  program  tdiich  implements  the  tedmique. 

Code's  method  first  uses  Carter's  results  for  products, 
quotients,  and  powers  so  that  a  sum  of  independent  H-function 
variates  remains.  The  Laplace  transform  of  each  term  in  the 
sun  is  eCtained,  then  evaluated  and  multiplied  at  corresponding 
values  of  the  transform  variable.  This  yields  a  tabular 
representation  for  the  Laplace  transform  of  the  sum.  This  is 
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then  numerically  inverted  iron  transform  space  using  Crunnp's 
method,  yielding  a  tabular  representation  of  the  p.d.f.  for  the 
algebraic  combination. 

Since  only  a  tabular  representation  of  the  p.d.f.  of  the 
sum  is  available  with  Codk's  method,  it  is  not  possible  to  use 
carter's  results  at  the  next  more  general  level  as  in  finding 
the  density  of  Z  in  Eg  (4.1).  Huts,  Cook's  method  will  work 
for  certain  algebraic  combinations  involving  products, 
quotients,  powers,  and  sums  of  ixKiependent  H-function  variates, 
but  not  others. 

Springer  [1979,  pp.  250-268]  describes  an  approach  due  to 
Carter  [1972]  to  approximate  the  p.d.f.  for  the  sum  or 
difference  of  independent  variates  based  on  the  moments  of  the 
sun  or  difference.  Carter  [1972]  also  wrote  a  FORTRAN  computer 
program  to  calculate  the  monents  of  an  algebraic  combination  of 
independent  H-function  variates  and  approximate  the  p.d.f.  and 
c.d.f.  from  these  moments.  Ibe  approximation  procedure  was 
developed  by  Hill  [1969]  and,  if  possible,  uses  either  a  Gram- 
Char  lier  type  A  series  (Hermite  polynomial)  or  a  Laguerre 
polynomial  series.  If  a  series  approximation  is  not  possible, 
the  first  four  moments  are  used  to  fit  a  probability 
distribution  from  the  Pearson  family.  As  Carter  [1972]  himself 
notes  "...  there  were  many  situations  in  tdiich  the  methods  did 
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not  work  or  in  vdiich  the  approxinations  were  totally 
unsatisfactory . " 

Since,  at  best,  only  a  series  approxination  of  the  p.d.f. 
of  the  sum  or  difference  is  available  with  this  method,  it  is 
still  not  possible  to  use  Carter's  results  at  the  next  more 
general  level  as  in  finding  the  density  of  Z  in  Eg  (4.1). 
'Hius,  this  method  will  also  work  for  certain  algdoraic 
conixLnations  involving  products,  quotients,  powers,  sums,  and 
differences  of  independent  H-function  variates,  but  not  others. 

4.1.  MOMBfTS  OF  TOE  SOM 

A  practical  method  to  approximate  the  density  of  the  sun 
of  independent  H-fxmction  variates  with  another  H-function 
distribution  is  based  on  the  moments  of  the  sum.  Ibe  moments 
of  the  sun  can  be  computed  using  the  moments  of  the  individual 
variates.  Svppose  X^,X2,...,X^  are  mutually  independent  random 
variables  and 
n 

Y  X.  (4.5) 

Ihe  r^^  moment  about  the  origin  of  Y  is 

=  (W  (<•«> 

For  the  case  n=2,  the  moments  of  Y  can  be  found  in  terms 
of  the  moments  of  and  X2  by  using  the  binomial  formula. 
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=  e[  ]  =  e[  (vEj)'  ] 


(4.7) 


since  mathenatical  expectation  is  a  linear  operator  and  and 
are  independent.  Here,  ^  f  j  is  the  bincndal  coefficient 
defined  as 


r! 

i!  (r-i)! 


(4.8) 


Eq  (4.7)  gives  the  mcments  of  T  in  terms  of  the  moments  of  X^ 
and  X^.  There  are  r+1  terms  in  the  sun  in  Eg  (4.7). 

For  the  case  0^2,  the  moments  of  Y  can  be  found  in  terms 
of  the  mcments  of  X^  to  X^  by  using  a  generalization  of  the 
binomial  formula. 


p_^  =  E[y^]  =e[  (vV 
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where  the  svm  is  over  all  i, such  that 

12  n 

n 

=  r  (4.10) 

Here,  the  multinonnial  coefficient  replaces  the  binomial 
coefficient.  Eq  (4.9)  gives  the  moments  of  Y  in  terms  of  the 
moments  of  through  There  are  ^  ^  j  terms  in  the  sum 

in  Eq  (4.9). 

4.2.  H-FTJNCTION  RARftMETER  ESTIMATES 

It  is  possible  to  use  the  moments  of  an  H-f unction  to 
estimate  its  parameters.  Bodenschatz  and  Boedigheimer  [1983; 
Boedigheimer  et  al,  1984]  devel(^>ed  and  verified  an  effective 
and  reliable  method  to  estimate  the  H-function  parameters  using 
the  method  of  moments.  The  tedmique  can  be  used  to  curve-fit 
a  mathematical  function  or  to  estimate  the  density  of  a 
particular  probability  distribution.  Tlieir  FORTRAN  computer 
program  will  accept  known  moments,  univariate  data,  ordered 
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pair  data  from  a  relative  frequency,  or  ordered  pair  data 
directly  from  the  function.  Output  from  the  program  are  the 
parameters  of  the  H-function  vdiose  moments  most  closely  match 
the  given  moments. 

This  technique  to  estinate  the  parameters  of  the  H- 
function  from  the  exact  moments  of  an  unknown  distribution  was 
also  generally  described  by  Jacdss  et  al  [1987].  A  more 
conplete  derivation  is  given  below,  tdiich  also  shows  how, 
through  algebraic  manipulation,  it  is  possible  to  reduce  the 
n\xiiber  of  equations  by  two. 

4.2.1.  METHOD  OF  MOKENTS 

The  method  of  moments  equates  an  appropriate  nurnber  of 
analytic  moments  of  the  H-function  with  corresponding  known 
moments  or  moments  calculated  from  data.  The  method  uses 
2(p>q)-*-2  moments  because  there  are  this  many  parameters  in  an 
H-ftmction  distribution.  This  produces  a  system  of  nonlinear 
equations  in  the  parameters  of  the  H-fxmction. 

Eq  (4.9)  gives  the  exact  moments  of  a  s\xn  of  independent 
variates  in  terms  of  the  moments  of  the  individual  variates. 
If  each  random  variable  in  the  sun  follows  an  H-function 
distribution,  the  exact  moments  of  the  individual  variates  are 
available  by  Eq  (3.5).  Using  these  moments  in  Eq  (4.9)  yields 
the  exact  moments  of  the  sun  of  independent  H-fuiction 
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variates . 

There  are  2(p-H3)+2  unknown  parameters  of  the  H-fxmction 
distribution,  specifically  (a^,A^),  i=l,...,p, 

j=l,...,q,  k,  and  c.  It  is  necessary,  therefore,  to  create 
2(p+q)+2  equations.  If  consecutive  moments  of  the  sun  are  vised 
starting  with  the  zeroth  mcment,  the  equations  are 

Mr  =  e[  ]  =  r^.1  I(r+1)  (4.11) 


for  r=0,l, . . . ,2(p+q)+l  where 


I(r+1)  = 


n  rfb.+B.+B.r) 
j=l  I  J  J  J  J 


n  r[a.+A.+A.r] 

i=n+l  ^  ^  ^  ^  J 


n 

n  r 

i=l 


q 

n  r 

j=mfl 


(l-bj-Bj-Bjr) 


(4.12) 

It  is  obvious  the  equations  are  nonlinear  since  they  involve 
gamma  functions.  Given  the  exact  mcments  of  the  svm,  the 
system  of  nonlinear  equations  as  in  Bq  (4.11)  needs  to  be 
solved  for  the  H-f unction  parameters. 


4.2.2.  REDUCING  TOE  SYSTai  OF  NCMLINEftR  EOGATICHS 
It  is  possible  to  eliminate  the  paurameters  k  and  c  from 
the  system  of  equations  through  algebraic  manipulation.  This 
procedure  reduces  by  two  the  number  of  simultaneous  nonlinear 
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equations  need  to  be  solved  to  give  the  H-functicn 

parameter  estimates.  Each  equaticn  of  Eq  (4.11)  can  be  solved 
for  k,  producing 

r+1 

Mr  c 

k  =  - - -  (4.13) 

I (r+1) 

for  r=0, . . .  ,2(p+q)+l,  where  I(r+1)  is  given  1^  Eq  (4.12). 
Since  all  the  equations  equal  k. 


k  = 


^0  g 
1(1) 


Ml  c 


1(2) 


Mj  c 
1(3) 


^2(p+q)+l 


2(p+q)+2 


l(2(p+q)+2) 


(4.14) 


tdiere  I(r+1)  is  given  by  Eq  (4.12).  The  adjacent  equations  in 
Eq  (4.14)  can  be  solved  for  c  to  give 

M^  I(r+2) 

c  =  - 1 -  (4.15) 

"t+i 


for  r=0, . . . ,2(p+q),  where  I(r+1)  is  given  by  Eq  (4.12).  Since 


all  the  equations  equal  c 
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/iQ  1(2)  Ml  1(3) 

Ml  1(1)  Mj  1(2) 


i^re  I(r+1)  is  given  by  Eg  (4.12).  The  adjacent  equations  in 
Eq  (4.16)  can  be  solved  to  give  the  following  honogeneous 
equations: 


^i  ^i+2 

(Mi+i)^  I(i+1)  X(i+3) 


(4.17) 


for  i=0,...,2(p+q)-l,  idxere  I(r+1)  is  given  by  Bq  (4.12).  This 
reduced  system  of  equations  involves  only  the  H-function 
parameters  (a ,A^)  ,  i=l,...,p.  and  ,  j=l,...,q. 

cnee  the  nonlinear  system  of  equations  in  Eq  (4.17)  is 
solved,  c  and  k  can  be  found  by  backsubstitution.  Using 
Bq  (4.16), 


c  = 

Ail  1(1) 


(4.18) 


vdiere  I(r+1)  is  given  by  Eq  (4.12).  Using  this  estimate  of  c 
in  Eq  (4.14)  produces 
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k=  — ? -  (4.19) 

1(1) 

vriiere  I(r-i-l)  is  givai  Eg  (4.12). 

Note  that  although  this  algebraic  nanipulation  has  reduced 
the  mxiber  of  equations  from  2(p+q)+2  to  2(p+q),  all  2(p+q)+2 
moments  of  the  sum  are  still  necessary. 

4.2.3.  SOLVING  TOE  SySTOi  OF  NQMLINERR  ECmTIOMS 

Even  after  reducing  the  system  of  nonlinear  equations  by 
two,  there  will  always  be  at  least  two  simultaneous  equations 
involving  gamma  functions.  To  complicate  matters  further,  the 
unknowns  (the  H-function  parameters)  appear  in  the  argunrent  of 
every  gamma  function.  A  general  analytic  solution  seems 
unlikely,  if  not  impossible. 

Bodenschatz  and  Boedigheimer  [1983]  conducted  a  literature 
review  of  various  nunerical  methods  applicable  to  the  problem. 
Most  methods  are  based  on  Newton's  method,  which  uses  the 
Jacobian  matrix  of  first  partial  derivatives  to  move  toward  the 
solution. 

Powell  developed  a  quasi -Newton  hybrid  algorithm  vdiich 
includes  the  beneficial  features  of  the  Levenberg-Marquardt 
method  and  implensnts  the  calculation-saving  strategy  of 
Broyden's  procedure.  Powell's  method  was  available  in  an 
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International  Mathematical  and  Statistical  Library  (IMSL) 
routine  named  ZSPGN.  In  conparison  to  available  software 
designed  to  solve  systems  of  nonlinear  equations,  ZSPOW  had 
outstanding  performance  and  initial  estimates  of  the  parameters 
had  little  effect  on  the  algorithm’s  stper  linear  convergence 
[Hiebert,  1980]. 

Bodenschatz  and  Boedigheimer  [1983]  used  ZSPGW  to  solve 
many  systems  of  nonlinear  equations  zis  in  Eq  (4.17)  with  great 
success.  In  their  FORTRAN  conputer  program,  they  allow  the 
user  to  make  initial  estimates  of  the  H-function  parameters  and 
also  provide  default  guesses  which  nearly  always  led  to 
convergence  to  the  correct  parameters. 

When  trying  to  fit  a  third-  or  higher-order  H-function, 
however,  the  nvinerical  solution  of  the  system  of  nonlinear 
equations  with  ZSPOW  is  sometimes  nunnerically  unstable. 
Several  other  numerical  methods  were  tried,  but  none  worked  as 
well  as  ZSPCW  in  solving  the  type  of  nonlinear  equations 
generated  using  the  method  of  moments  with  the  H-function.  An 
area  of  further  research  is  to  develop  a  better  way  to  solve 
the  system  of  nonlinear  equations. 

Jacobs  et  al  [1987]  discussed  using  equations  as  in 
Eq  (4.11)  to  estimate  the  parameters  of  an  H-function  if  the 
exact  moments,  of  the  unknown  distribution  were  available. 
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They  also  gave  analytic  solutions  to  these  equations  for 
certain  first  order  H-function  distributions  with  B=1  or 
B=  — However,  restrictions  such  as  these  nay  not  give  an  H- 
function  with  as  good  a  fit  to  the  moments  eis  an  unrestricted 
solution  to  the  equations  might  allow. 

*4.3.  SPECIAL  OC^SIDERATIOIS  FCSt  TYPE  VI  H-FONCTICM  VaRIATES 

Suns  of  certain  H-function  distributions  have  special 

properties  which  allow  an  analytic  solution  to  the  system  of 

n 

nonlinear  equations.  For  exarrple,  if  Y  =  V  X.  vdiere  the  X. 

i=l  ^  ^ 

are  nutiially  independent  random  variables  which  have  beta 

distributions  over  (0,1),  then  the  p.d.f.  of  Y  will  be  nonzero 

n 

only  over  (0,n).  Similarly,  if  Y  =  T  X.  where  the  X.  are 

i=l 

mutually  independent  random  variables  which  have  Pareto 
distributions  over  (l,<x)),  then  the  p.d.f.  of  Y  will  be  nonzero 
only  over  (n,a>).  A  new  way  to  take  advcuitage  of  these  ranges 
was  discovered  which  eliminates  the  need  for  a  nunerical 
solution  to  the  system  of  nonlinear  equations. 

Certain  siin)lifications  occur  if  each  H-fmction 
distribution  in  the  s\in  is  a  Type  VI  convergent  H-function 
according  to  Table  1  in  Section  2.3.  Because  the  suns 
described  above  will  also  have  restricted  ranges,  they  must 
also  be  Type  VI  convergent  H-f unctions.  Sections  2.6  and  3.6 
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gave  the  new  H-f\mction  representatioons  of  these  types  of 
f  met  ions  over  more  general  ranges. 

J  Type  VI  convergent  H-fmetions  have  a  "B"  ganma  term 
in  the  nimerator  and  an  "A"  gamra  term  in  the  denominator  of 
the  definition  of  the  H-fmcticn,  with  A=B.  Many  ^  Type  VI 
convergent  H-f motions  have  A=B=1.  Scaling  over  a  range 
different  from  (0,1)  is  achieved  throxigh  the  H-fmetion 
parameter  c.  If  the  range  of  the  srm,  Y,  is  only  over  (0,n), 
then  in  the  H-f  motion  representation  of  the  p.d.f.  of  Y, 
c=  Exploiting  the  imiqiie  nature  of  ^  Type  VI 
convergent  H-fmetions  has  reduced  the  nurber  of  parameters  to 
estimate  from  six  to  three.  The  three  nonlinear  equations  in 
the  ixnknown  parameters  k,  b,  and  a  are 


k 

r(b+i) 

c 

r(a+l) 

k 

r(b+2) 

c^ 

r(a+2) 

k 

r(b+3) 

c^ 

r(a+3) 

(4.20) 


(4.21) 


(4.22) 


where  c=  the  range  of  the  sum  is  over  (0,n),  and 

through  ^2  exact  moments  of  the  sun  from  Eq  4.9.  These 

three  nonlinear  equations  allow  an  analytic  solution.  The 


solution  is 
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a  = 


^0 


^0  *^2 


[V2  -  (»^l)  ^0^2  ■  t'^l] 


-  1 


* 


* 


b  = 


k  = 


Mq  c 


r(a+i) 

r(b+i) 


-  1 


(4.23) 

(4.24) 

(4.25) 


where  c  and  throiigh  are  as  given  above.  A  possible  H- 
fimction  representation  for  the  density  of  the  stm  of  ^  Type 
VI  convergent  H-functions  is 

*  f(y)  =  k  ^  y  :  (a,i)  ;  (b,l)  ]  (4.26) 

where  the  range  of  the  sum  is  over  (0,n)  and  a,  b,  and  k  are 
given  by  Eq  (4.23),  Eq  (4.24),  and  Eq  (4.25),  respectively. 

A  similar  approach  applies  to  sixns  of  ^  Type  VI 
convergent  H-functions  which  have  an  "A”  ganna  term  in  the 
nvmerator  and  a  "B"  gamta  term  in  the  denominator  of  the 
definition  of  the  H-f unction,  with  A=B.  Again,  let  A=B=1. 
Scaling  over  a  range  different  from  (l,oo)  is  achieved  through 
the  H-f  unction  parameter  c.  If  the  range  of  the  sum,  Y,  is 
over  (n,oD),  then  in  the  H-fvinction  representation  of  the  p.d.f. 
of  Y,  c=  Exploiting  the  unique  nature  of  H^  ^  Type  VI 

convergent  H-functions,  the  three  nonlinear  equations  in  the 
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unknown  parameters  k,  b,  and  a  are 

k 

Mft  =  -rr- 

r(-a) 

(4.27) 

0  c 

r(-b) 

k 

/i-i  —  — ij— 

r(-a-l) 

(4.28) 

1  2 
c 

r(-b-i) 

k 

^2  ■  — 3“ 

r(-a-2) 

(4.29) 

r(-b-2) 

where  c=  the  range  of  the  stm  is  over  (n,oc>),  and  /iQ 

throijgh  are  the  exact  moments  of  the  stm  from  Eq  4.9.  These 
three  nonlinear  equations  again  allow  an  analytic  solution. 
The  solution  for  the  unknown  parameters  a,  b,  and  k  is 
identical  to  the  ones  given  by  Eq  (4.23),  Eq  (4.24),  ajnd 
Eq  (4.25),  respectively.  A  possible  H-function  representation 
for  the  density  of  the  sum  of  ^  Type  VI  convergent  H- 
f unctions  is 

*  f(y)  =  k  hJ  ^  y  :  (a,l)  ;  (b,l)  ]  (4.30) 

where  the  range  of  the  sun  is  over  (n,oo)  and  a,  b,  and  k  are 
given  by  Eq  (4.23),  Eq  (4.24),  and  Eq  (4.25),  respectively. 

*4.4.  DElOtSTRATION  OF  THE  TECHNIQUE 

A  FCXTTRAN  computer  program  was  developed  to  inplement  the 
new  technique  of  finding  an  H-function  distribution  which 
approximates  the  distribution  of  the  sun  of  two  or  more 
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indepeafident  H-f unction  variates.  The  ccnnputer  program  was 
designed  to  run  interactively,  asking  questions  of  the  user  and 
expecting  a  response. 

Input  to  the  program  can  be  done  interactively  or  through 
an  input  file.  A  user  will  input  the  n\xii)er  of  terms 
(independent  random  variables)  in  the  sun  and  the  type  of  each 
variate.  As  currently  configured,  the  program  will  accept  V5> 
to  five  terms  in  the  sum.  The  type  of  each  variate  may  be  one 
of  the  named  special  cases  or  a  general  H-function 
distribution.  If  the  type  is  a  special  case,  the  program  will 
convert  the  parameters  of  the  special  case  into  the  correct  H- 
f  unction  parameters  for  the  H-function  r^resentation. 
Alternatively,  a  \iser  may  input  the  H-function  parameters 
directly. 

The  program  will  then  conpute,  using  Bq  (3.5),  the  correct 
moments  of  each  H-fvnction  variate  in  the  sun.  The  program 
will  verify  a  zeroth  moment  of  unity  and  query  the  \iser  if  the 
H-function  parameters  do  not  give  a  valid  density.  If  desired, 
the  program  will  reconpute  the  constant  k  to  produce  a  valid  H- 
f unction  distribution.  The  moments  of  each  term  in  the  sun  are 
used  as  in  Eq  (4.9)  to  conpute  the  exact  moments  of  the  sun. 

If  the  sum  involves  only  H^  ^  Type  VI  H-function 
distributions  or  H^  ^  Type  VI  H-function  distributions,  the 
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analytic  solutions  for  the  H-function  parameters  derived  in 
Section  4.3  are  used.  Otherwise,  an  adapted  version  of  the 
program  of  Bodenschatz  and  Boedigheimer  [1983]  is  called  to 
create  and  solve,  using  ZSPOW,  the  system  of  nonlinear 
equations  as  in  Eq  (4.17).  This  produces  parameters  of  the  H- 
f unction  distribution  whose  moments  will  closely  match  the 
exact  moments  of  the  sum. 

If  there  are  only  two  terms  in  the  sun,  an  attempt  can  be 
made  to  fit  an  H-function  distribution  with  tp  to  five  gartma 
terms  (fifth-order).  If  there  are  more  than  two  terms  in  the 
sum,  the  program  currently  computes  only  the  zeroth  through  the 
fifth  moments  about  the  origin  of  the  sum.  This  limits  the 
choice  of  an  H-ftmction  distribution  for  the  sum  to  a  first-  or 
second-order  H-function. 

This  is  not  a  serious  restriction  for  two  reasons.  First- 
and  second-order  H-function  distributions  can  represent  a  wide 
variety  of  distributional  forms,  including  nearly  every  named 
special  case.  Secondly,  nunerical  eund  computational 
lindtations  exist  in  finding  the  moments  of  the  sum  and  solving 
large  systems  of  nonlinear  equations  with  ZSFOH.  Until  these 
limitations  are  addressed,  an  attempt  to  fit  a  third-  or 
higher-order  H-function  would  probably  not  be  successful . 

The  moments  of  the  new  H-function  distribution  are 
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canrputed  for  conparison  to  the  exact  moments  of  the  sum.  The 
program  will  then,  if  desired,  create  an  input  file  for  another 
conputer  program  [Ck)ok,  1981;  Code  and  Barnes,  1981]  which  will 
evaluate  and  plot  the  resulting  p.d.f.  and  c.d.f. 

An  output  file  is  always  created,  giving  the  H-functicn 
parameters  and  moments  of  each  term  in  the  sixn,  the  exact 
moments  of  the  sxxn,  the  estimated  H-f motion  parameters,  and 
the  moments  of  this  new  H-f unction  distribution. 

The  PCaWRAN  source  code  of  the  corputer  program  described 
eJsove  is  not  provided  in  an  appendix.  Instead,  it  is  avadlable 
from  Dr.  J.  Wesley  Barnes  in  the  Department  of  Mechanical 
Engineering  at  the  Uhiversity  of  Texas  at  Austin.  If 
requested.  Dr.  Barnes  will  transfer  the  program  by  electronic 
mail  or  by  floppy  disk. 

The  new  FORTRAN  conputer  program  which  iirplements  the 
technique  currently  operates  on  the  instructional  VAX  cluster 
at  the  University  of  Texas  at  Austin.  This  cluster  links  two 
VAX  6420  cenputers  and  several  VAX  11/780  conputers  and  uses 
the  VMS  operating  system.  The  VAX  PCSWRAN  cotipiler  was  used 
and  the  program  linked  to  IMSL  for  the  call  to  ZSPOW.  The 
program  should  successfully  run  on  other  computers  which 
support  ZSPOW  in  IMSL  and  which  have  ANSI  FORTRAN  capabilities. 

The  new  technique  and  conputer  program  are  demonstrated  in 
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the  following  exanples.  In  all  but  the  final  example,  it  is 
possible  to  find  the  exact  distribution  of  these  sums  of 
variates  by  other  methods.  Thus,  we  can  measure  the 
effectiveness  of  the  new  technique  by  comparing  the  resulting 
H- function  distribution  to  the  exact  distribution.  Of  course, 
the  new  technique  also  works  with  sixns  for  which  the  solution 
by  other  methods  is  very  difficult  (or  ini>ossible  for  all 
practical  purposes). 

4.4.1.  EXaMPLE  1  -  SUM  OF  TOREE  INDEPEMDENT.  IDElfriCALLY 

DISTRIBUTED  GAMMA  VARIATES 

Suppose  we  want  to  use  this  tedmique  to  find  the  H- 
f unction  distribution  of  the  sum  of  three  independent, 
identically  distributed  ganma  variates  with  parameters  r=2  and 
X=2.  Using  the  computer  program,  we  input  ”3"  as  the  nuniber  of 
terms  in  the  sum,  "ganrna"  as  the  type  of  each  variate,  and  the 
parameters  "2"  for  r  and  "2”  for  X. 

The  momients  of  each  variate  in  the  sum  are  computed  and 
combined  to  give  the  exact  moments  of  the  sixn.  These  are 
/Jq  =  1.0000000  Ml  =  3.0000000 

Mj  =  10.5000000  =  42.0000000 

M4  =  189.0000000  M5  =  945.0000000 

The  program  asks  whether  we  want  to  fit  an  H-fvjnction 
distribution  to  the  moments  of  the  s\jm  and,  if  so,  the  order  of 
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H-function  to  fit.  Since  all  three  variates  in  the  siin  are 
first  order  H-functions,  we  try  to  fit  another  first  order  H- 
f unction.  We  enter  "1  0  0  1"  for  m,  n,  p,  and  q,  respectively. 
The  program  uses  the  zeroth  throxigh  third  moments  of  the  sum  to 
create  and  solve  the  system  of  nonlinear  equations,  producing 
the  H-function  parameter  estimates 
b  =  4.99999988853 
B  =  0.99999999256 
k  =  0.01666666667 
c  =  2.00000000000 

Allowing  for  roundoff  error,  we  recognize  these  ais  the 
parameters  of  the  H-function  representation  of  the  gamma 
distribution  with  r=6  and  X=2.  The  procedure  found  the 
"correct"  distribution  since  the  gamma  distribution  with  a 
conmon  X  has  the  reproductive  property.  The  program 
automatically  uses  these  H-function  parameters  to  ccnpute  the 
moments  of  the  new  H-function  distribution. 

/Jq  =  1.0000000  =  3.0000000 

/i2  =  10.5000000  ^3  =  42.0000000 

=  189.0000000  ^5  =  945.0000000 

The  program  then  asks  whether  we  want  to  create  an  input 
file  for  another  cotrputer  program  [Cook,  1981;  Cook  and  Barnes, 
1981]  which  will  evaluate  and  plot  the  resulting  p.d.f.  and 
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c.d.f.  If  so,  it  oisks  the  range  and  interval  of  the  variable 
where  we  want  these  fmctions  eval\iated.  We  try  the  range  0.1 
to  7.0  with  an  interval  of  0.1.  The  ccrplete  output  file  is 
provided  in  Appendix  A. 


4.4.2.  EXAMPLE  2  -  SUM  OF  TWO  INDEPEMDEUT  ERLANG  VARIATES 


WITH  DIFb'EkENT  X 


Suppose  we  want  to  find  an  H-fmction  distribution  for 
Y  =  ^^+^2  where  has  an  Erlang  (gaima)  distribution  with 
parameters  r=2  and  X=4  and  X2  has  an  Erlang  (ganna) 
distribution  with  parameters  r=l  and  X=2.  Note  that  X2  also 
has  an  exponential  distribution  with  X=2.  Using  the  conputer 
program,  we  input  "2"  as  the  number  of  terms  in  the  sum, 
"ganma"  as  the  type  of  each  variate,  the  parameters  r="2"  and 
X="4"  for  X^,  and  the  parameters  r=”l'’  and  X="2"  for  X2.  Of 
course,  we  could  have  input  X2  cus  type  "exponential"  with  X="2" 
or  input  both  random  variables  as  general  (unnamed)  H-function 
variates  by  giving  their  H-function  parameter  representations. 

The  moments  of  each  variate  in  the  sum  are  computed  and 
conribined  to  give  the  exact  moments  of  the  sum.  These  are 


Mg  =  1.0000000 

=  1.3750000 
=  5.3437500 

Mg  =  43.4179688 


=  1.0000000 
M3  =  2.4375000 

Mg  =  14.0625000 

M7  =  154.4238281 
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Since  both  variates  in  the  sum  are  first  order  H- 
f unctions,  we  try  to  fit  another  first  order  H-fmction  to  the 
moments  of  the  sum.  We  enter  "1  0  0  1"  for  m,  n,  p,  and  q, 
respectively.  The  program  uses  the  zeroth  through  third 
moments  of  the  sum  to  create  and  solve  the  system  of  nonlinear 
equations,  producing  the  H-f unction  parameter  estimates 
b  =  3.02485587519 
B  =  1.26389176312 
k  =  0.74895896608 
c  =  6.53330240202 

The  program  automatically  uses  these  H-function  parameters 
to  compute  the  moments  of  the  new  H-function  distribution. 

Mq  =  1.0000000  =  1.0000000 

=  1.3750000  /J3  =  2.4375000 

=  5.3382458  =  14.0112602 

Mg  =  43.0835662  M7  =  152.4885870 

We  choose  to  evaluate  and  plot  the  resulting  p.d.f.  and 
c.d.f.  over  the  range  0.05  to  5.50  with  an  interval  of  0.05. 
The  complete  output  file  is  provided  in  Appendix  A. 

4.4.3.  EXAMPLE  3  -  SUM  OF  TWO  INDEPEMDEUT  STANDARD 

UWIFCM4  VARIATES 

Suppose  we  want  to  find  an  H-function  distribution  for 
Y  =  where  and  X2  each  have  a  uniform  distribution 
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over  (0,1).  Using  the  carputer  program,  we  input  "2"  as  the 
number  of  terms  in  the  sum  and  ’’taiiform  over  (0,1)"  as  the  type 
of  each  variate. 

The  moments  of  each  variate  in  the  sun  are  confuted  and 
cciTi>ined  to  give  the  exact  moments  of  the  sun.  These  are 


/Jq  =  1.0000000 

=  1.1666667 

=  2.0666667 

Pg  =  4.5357143 

Pg  =  11.3555556 

=  31.0151515 


=  1.0000000 
/ig  =  1.5000000 

/jg  =  3.0000000 

=  7.0833333 

Mg  =  18.6000000 


Since  both  variates  in  the  sum  are  Type  VI  convergent  J 
H- functions  over  (0,1),  the  program  uses  the  zeroth  through 
second  moments  of  the  sun  to  ccnpute  the  analytic  solution  for 
the  H-f unction  parameters.  We  get  the  H-function  paraometer 
estimates 


b  =  1.50000000000 
B  =  1.00000000000 
a  =  4.00000000000 
A  =  1.00000000000 
k  =  9.02703333685 


c  =  0.50000000000 

The  program  autcmatically  uses  these  H-function  parameters 
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to  connpute  the  monents  of  the  new  H-f unction  distribution. 


*^0  = 

1.0000000 

=  1.0000000 

1.1666667 

^2  =  1.5000000 

^4  = 

2.0625000 

=  2.9791667 

= 

4.4687500 

=  6.9062500 

II 

00 

a. 

10.9348958 

Mg  =  17.6640625 

*^10  = 

29.0195312 

We  choose  to  evaluate  and  plot  the  resulting  p.d.f.  and 
c.d.f.  over  the  range  0.05  to  2.00  with  an  interval  of  0.05. 
The  carpi ete  output  file  is  provided  in  Appendix  A. 

4.4.4.  EXAMPLE  4  -  SUM  OF  TWO  INDEPEM>EJrr .  IDEHTICALLY 

DISTRIBUTED  BETA  VARIATES 

Suppose  we  want  to  find  an  H-function  distribution  for 
Y  =  where  and  X2  are  independent  variates  with  beta 

distributions  over  (0,1)  with  parameters  a=l  and  p=2.  Using 
the  corputer  program,  we  input  ”2"  as  the  nunriber  of  terms  in 
the  sum,  "beta  over  (0,1)"  as  the  type  of  each  variate,  "1"  for 
a,  and  "2"  for  p. 

The  moments  of  each  variate  in  the  son  are  conputed  and 
combined  to  give  the  exact  moments  of  the  sum.  These  are 
=  1.0000000  =  0.6666667 

^2  =  0.5555556  =  0.5333333 

^4  = 


0.5666667 


=  0.6507937 
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Mg  =  0.7952381  M7  =  1.0222222 
Mg  =  1.3703704  Mg  =  1.9030303 
Mj^q  =  2.7229437 

Since  both  variates  in  the  sum  are  Type  VI  convergent  ^ 
H- functions  over  (0,1),  the  program  uses  the  zeroth  through 
second  moments  of  the  sixn  to  compute  the  analytic  solution  for 
the  H-f unction  parameters.  We  get  the  H-f unction  parameter 
estimates 

b  =  1 . 33333333333 
B  =  1.00000000000 
a  =  6.00000000000 
A  =  1.00000000000 
k  =  302.35856086513 
c  =  0.50000000000 

The  program  automatically  uses  these  H-f  unction  parameters 
to  compute  the  moments  of  the  new  H-function  distribution. 


*^0  = 

1.0000000 

=  0.6666667 

0.5555556 

M3  =  0.5349794 

^4  = 

0.5706447 

Mg  =  0.6571060 

= 

0.8031296 

M7  =  1.0296533 

M 

00 

a. 

1.3728711 

Mg  =  1.8915113 

Mm  = 

2.6796410 

We  choose  to  evaluate  and  plot  the  resulting  p.d.f.  and 
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c.d.f.  over  the  range  0.05  to  2.00  with  an  interval  of  0.05. 
The  connplete  output  file  is  provided  in  Appendix  A. 


4.4.5.  EXAMPLE  5  -  SUM  OF  TWO  INDEPENDENT  VARIATES  WITH 


WEIBULL  AND  RAYLEIGH  DISTRIBUTICTIS 


A  final  exairple  demonstrates  the  technique  and  corputer 
program  when  the  exact  distribution  of  the  stnn  is  very 
difficult  or  impossible  to  obtain. 

Svppose  we  want  to  find  an  H-function  distribution  for 
Y  =  ^^+^2  where  has  a  Weibull  distribution  with  parameters 
^=5  and  X=4  and  X2  has  a  Rayleigh  distribution  with  parameter 
X=3.  Using  the  computer  program,  we  input  "2"  as  the  nxanber  of 
terms  in  the  sum,  "Weibull”  as  the  type  of  variate  X^^  with 
parameters  p="5"  and  X="4",  and  "Rayleigh"  eis  the  type  of 
variate  X2  with  parameter  X="3".  Of  course,  we  could  have 
input  both  random  variables  as  general  H-function  variates  by 
giving  their  H-function  parameter  representations. 

The  moments  of  each  variate  in  the  sum  are  computed  and 
combined  to  give  the  exact  moments  of  the  sum.  These  are 


/Jq  =  1.0000000 

^2  =  1.5550061 

=  3.0567312 

Mg  =  7.3134235 

Mg  =  20.6900000 


=  1.2075051 
M3  =  2.1228284 
Mg  =  4.6225066 
M7  =  12.0647384 
Mg  =  36.7858797 
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=  67.6453356 

Since  both  variates  in  the  stJtn  are  first  order  H- 
ftmctions,  we  try  to  fit  another  first  order  H-f unction  to  the 
moments  of  the  sum.  We  enter  "1  0  0  1"  for  m,  n,  p,  and  q, 
respectively.  The  program  uses  the  zeroth  through  third 
moments  of  the  sum  to  create  and  solve  the  system  of  nonlinear 
equations,  producing  the  H-f unction  parameter  estimates 
b  =  5.29104742893 
B  =  0.62358286635 
k  =  0.02369785777 
c  =  2.45993399245 

The  program  automatically  uses  these  H-f  unction  parameters 
to  compute  the  moments  of  the  new  H-f  unction  distribution. 


II 

o 

3. 

1.0000000 

=  1.2075051 

^2  = 

1.5550061 

^3  =  2.1228284 

*^4  = 

3.0568212 

/J3  =  4.6235079 

*^6  = 

7.3191453 

=  12.0889706 

II 

00 

3. 

20.7770136 

Mg  =  37.0686930 

*^10  = 

68.5084603 

We  choose  to  evaluate  and  plot  the  resulting  p.d.f.  and 
c.d.f.  over  the  range  0.05  to  2.40  with  an  interval  of  0.05. 
The  complete  output  file  is  provided  in  Appendix  A. 


CHAPTER  5 


CXMPARING  THE  ESTIMATED  H-roNCTIC»J  TO  THE 
EXACT  DISTRIBUTIOI  OF  THE  SIW 

It  is  natiiral  to  want  to  measure  the  effectiveness  of  the 
newly  developed  technique  which  finds  an  H-f unction 
distribution  whose  moments  closely  approximate  the  exact 
moments  of  the  sum  of  independent  H-f  unction  variates. 
Ideally,  the  approximate  H-f unction  distribution  should  be 
compared  to  the  exact  distribution  for  the  sim.  Unforttjnately, 
the  exact  distribution  of  the  sum  is  often  difficult  to  obtain. 
Because  of  this,  analysts  usually  resort  to  computer  simulation 
to  arialyze  the  resulting  distribution.  Even  this  approach 
yields  only  information  about  the  resulting  distribution,  not 
the  exact  density  itself. 

Four  examples  included  in  this  thesis  were  chosen  as 
representative  of  sunns  of  independent  H-function  variates  where 
the  t>xact  distribution  of  the  sum  can  be  obtained.  A  FQRTHAN 
computer  program  was  written  to  compare  the  approximate  H- 
function  distribution  to  the  exact  distribution  for  the  sxin, 
using  several  different  measures  of  merit.  The  newly  developed 
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technique  was  also  effective  in  finding  an  H-function 
distribution  whose  moments  closely  approximate  the  exact 
moments  of  the  sum  even  when  the  exact  distribution  is  very 
difficult  or  impossible  to  obtain. 

5.1.  FINDING  THE  EXACT  DISTRIBOTICW  OF  THE  SUM 

As  mentioned  above,  it  is  often  difficult  to  obtain  the 
exact  distribution  for  the  sxm  of  independent  random  variables. 
While  the  moment  generating  function,  Laplace  transform,  or 
characteristic  function  (Fourier  transform)  of  specific 
densities  are  sometimes  available  from  tables,  the  product  of 
these  transform  functions  will  often  produce  a  functional  form 
that  is  not  easily  inverted. 

Si9>pose  Xj^,X2, .  • .  ,X^  are  mutually  independent  random 
variables  with  respective  densities  fj^(Xj^) ,  f2(x2) , 

fn(xj  ,  each  nonzero  only  for  positive  values  of  the  varizd>le. 

n 

Suppose  we  want  the  density  of  Y  =  X.  .  If  the  Laplace 

i=l  ^ 

transforms  of  the  densities  are  Zg^f^^ 

^sKC^'n)}  ,  respectively,  then  the  density  of  Y  is  available  as 


the  inversion  integral 
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u+ioo 

p 

n 

r 

• 

^  ■SFI 

1 

*s{  ^C-i)  } 

ds 

(5.1) 

'^u-i® 

i=l 

• 

■ 

The  contour  integral  of  a  conplex  variable  in  Eq  (5.1)  is 
usually  difficult  to  evaluate.  TTie  approach  is  similar  (and  as 
difficult)  if  the  moment  generating  function  or  characteristic 
function  (Fourier  transform)  is  used  instead  of  the  Laplace 
transform. 

In  certain  situations,  there  are  ezisier  ways  to  find  the 
exact  distribution  of  the  svm  of  independent  random  variables. 
These  are  discussed  in  the  following  sections. 

5.1.1.  CXINVOLUTIOW  INTEX3RAL 

If  the  product  of  two  transform  functions  can  not  be 
recognized  as  the  transform  of  a  specific  function,  the 
inversion  may  be  done  with  the  convolution  integral .  If  and 
X2  are  independent  random  variables  with  respective  densities 
fl(xi)  and  £2(^2) '  nonzero  only  for  positive  values  of  the 
variable,  then  the  density  of  Y=Xj^+X2  is 

The  convolution  integral  in  Eq  (5.2)  is  often  eaisier  to 
evaluate  than  Eq  (5.1). 

The  convolution  integral  was  used  to  find  the  exact 
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distribution  of  the  sim  of  selected  random  variables  with 
uniform,  power  fimction,  or  beta  distributions.  Only  2  of  the 
15  exact  distributions  given  below  for  the  sum  of  two 
independent  variates  are  used  as  exarnples  in  this  thesis.  The 
others  are  only  included  for  the  reader's  benefit,  perhaps  by 
saving  the  reader  the  work  of  deriving  them.  Consider  the 
mutually  indeper-dent  random  variables  X^,X2, . . .  with 
densities 


0  <  Xj^  <  1 

(5.3) 

otherwise 


0  <  X2  <  1 

otherwise 


(5.4) 


0  <  Xj  <  1 

otherwise 

0  <  X.  <  1 

4 

otherwise 
0  <  x^  <  1 


(5.5) 


(5.6) 


^5W  = 


0 


otherwise 


(5.7) 
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‘7W  = 


2  -  2x, 


fgCXg)  = 


3|x 


^ioC*io)  ■  1 


f  5(i-»io) 


0  <  Xj  <  1 


otherwise 

0  <  x^  <  1 
otherwise 


0  <  Xg  <  1 


otherwise 


0  <  Xg  <  1 


otherwise 


0  <  'lo  <  1 


otherwise 


If  y=Xj^+Xg,  then  the  exact  p.d.f.  of  Y  is 

y  0  <  y  i  1 

f(y)  =  j  2-y  1  s  y  <  2 

0  otherwise 


(5.8) 


(5.9) 


(5.10) 


(5.11) 


(5.12) 


(5.13) 
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If  Y=X^+X2,  then  the  exact  p.d.f.  of  Y  is 


f(y)  = 


2y-y" 

0 


0  <  y  i  1 
1  i  y  <  2 
otherwise 


If  Y=X^+X^,  then  the  exact  p.d.f.  of  Y  is 
2 


f(y)  = 


r  2y-y" 
(y-2)‘ 
0 


0  <  y  ^  1 
1  s  y  <  2 
otherwise 


If  Y-Y.2^^1  then  the  exact  p.d.f.  of  Y  is 

0  <  y  s:  1 


f(y)  = 


2  y3 

-T  y 


■2  3  ^  .  8 

-r-  Y  +  4y  - 


1  s  y  <  2 
otherwise 


If  y=X2+Xg,  then  the  exact  p.d.f.  of  Y  is 


f(y)  =  \ 


y3  -  4y2  +  4y 

(y-2)^ 


0  <  y  s  1 
1  s  y  <  2 


(5.14) 


(5.15) 


(5.16) 


(5.17) 


0 


otherwise 
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Y=X2+X2,  then  the  exact  p.d.f.  of  Y  is 


f(y)  = 


^  2y^ 


4  y3  -  2/  + 


0  <  y  i  1 
1  i  y  <  2 
otherwise 


(5.18) 


If  Y=Xj^+X^,  then  the  exact  p.d.f.  of  Y  is 


f(y)  =  \  1  -  (y-i)' 


0  <  y  i  1 
1  i  y  <  2 
otherwise 


(5.19) 


If  Y=Xj^+Xg,  then  the  exact  p.d.f.  of  Y  is 


’  1  +  (y-1)' 

f(y)  =  .  -(y-2)^ 


0  <  y  i  1 
1  i  y  <  2 
otherwise 


(5.20) 


If  Y=X2+X^,  then  the  exact  p.d.f.  of  Y  is 


f(y)  = 


1  4 

~r  y 

y^  +  3y^  -  2y 


0  <  y  i  1 
1  s  y  <  2 
otherwise 


(5.21) 


If  y=X2+Xg,  then  the  exact  p.d.f.  of  Y  is 


f(y)  = 


-  2y^  +  3y^ 


0  <  y  i  1 


iy^  +  2y^  -  8y  +  8  1  i  y  <  2 


otherwise 


(5.22) 


If  Y=X2+X^,  then  the  exact  p.d.f.  of  Y  is 


f(y)  = 


-1  4^.3 

—>r-  y  +  2y 


14  3  2 

y  “  ^y^  -  4y  +  4 


0  <  y  i  1 

1  i  y  <  2  (5.23) 

otherwise 


If  y=X2+Xg,  then  the  exact  p.d.f.  of  Y  is 


-14  3  2 

-4-  y  +  4y'’  -  9y^  +  6y 


f(y)  =  J  4.  (y-2)4 


0  <  y  i  1 

1  s  y  <  2  (5.24) 

otherwise 


If  Y=X^+Xg,  then  the  exact  p.d.f.  of  Y  is 

,3  5 

TiT  ^ 

£(y)  =  ■  y®  +  6y2  -  9y  + 


0  <  y  i  1 
1  s  y  <  2  (5.25) 


otherwise 
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If  y=Xc+X,rt/  then  the  exact  p.d.f.  of  Y  is 

D  lU 


f(y)  = 


3 

“10" 

-3 

ID" 


•i  4  2 

y  -  3y  +  12y  -  18y^  +  9y 

(y-2)^ 


0  <  y  i  1 
1  i  y  <  2 
otherwise 


(5.26) 


If  Y=X.+X_,  then  the  exact  p.d.f.  of  Y  is 
4  D 


f(y)  = 


3  4^,3 

y  +  3y 


3  5 

-^y^  +  -|-y^  -  3y^  +  6y2  -  12y  + 


48 


0  <  y  i  1 
1  i  y  <  2 
otherwise 
(5.27) 


These  exact  distributions  could  be  cotrpared  to  the 
corresponding  approximate  H-f unction  distributions.  The  exact 
distribution  in  Eq  (5.13)  is  connpared  to  the  approximate  H- 
function  distribution  in  Section  5.3.3  below.  The  exact 
distribution  in  Eq  (5.17)  is  corpared  to  the  approximate  H- 
function  distribution  in  Section  5.3.4  below. 

Appendix  B  contains  graphical  depictions  of  the  15  exact 
distributions  in  Eq  (5.13)  through  (5.27)  zJbove.  The  graph  of 
the  p.d.f.  of  each  random  varicdsle  in  the  sum  is  followed  by 
the  graph  of  the  p.d.f.  of  the  sum.  Each  series  of  three 
graphs  is  referenced  to  the  corresponding  equation  nxxrber. 
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5.1.2.  REPRODUCTIVE  DISTRIBOTICWS 

As  discussed  in  Secticaa  1.4. 1.3,  certain  distributions 
have  the  reproductive  property,  which  considerably  simplifies 
finding  the  exact  distribution  of  the  sum  of  independent  random 
variables  with  distributions  of  these  forms.  The  ganma 
distribution  is  an  H-function  distribution  with  the 
reproductive  property,  provided  X  is  common  among  all 
distributions  in  the  sum. 

n 

Suppose  we  want  the  density  of  Y  =  X.  where 

i=l  ^ 

X^,X2,...,X^  are  mutually  independent  gamma  distributed  random 
variables  with  respective  parameters  (rj^,X) ,  (r2A) »  ...» 

(rn»X)  ,  and  X  is  common  to  all  the  distributions.  It  is  well 
known  that  Y  also  has  a  gamnna  distribution  with  parameters 
n 

r  =  ^  r.  and  X.  Since  the  exponential  distribution  is  a 
i=l  ^ 

special  case  of  the  ganma  distribution  with  r=l,  the  svm  of  n 
independent,  identically  distributed  exponential  random 
variables  has  a  gamma  distribution  with  parameters  r=n  and  X. 
This  result  is  also  widely  known. 

The  reproductive  property  of  the  gamma  distribution  was 
used  to  find  the  exact  distribution  of  the  sum  of  three 
independent,  identically  distributed  gamma  variates.  This 
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exact  distribution  is  compared  to  the  approximate  H-function 
distribution  in  Section  5.3.1  below. 


*  5.1.3.  ERLANG  DISTRIBOTIC»>S  WITH  DIFFERBn'  X 

Another  new  finding  permits  finding  the  exact  distribution 
of  the  sum  of  independent  Erlang  distributed  random  variables, 
even  when  X  is  different.  An  Erlang  distribution  is  simply  a 
ganrna  distribution  with  an  integer  shape  parameter  r. 

Suppose  we  want  the  density  of  y=X^+X2  where  X^  and  X2  are 
independent  Erlang  distributed  random  varicd>les  with  respective 
parameters  (r^,X^)  and  (r2/X2)  /  with  Xj^^X2.  The  Laplace 
transform  of  the  density  of  Y  is  the  product  of  the  Laplace 
transforms  of  the  Erlang  densities. 


(5.23) 


We  proceed  by  decomposing  the  final  two  terms  by  partial 
fractions.  The  approach  presented  below  is  essentially  the 
same  as  that  described  by  Kleinrock  [1975]  and  is  routinely 
given  in  texts  covering  Laplace  transforms  inclviding  Churchill 
[1972],  Widder  [1941],  Thompson  [1960],  Shiith  [1966],  Doetsch 
[1971],  Davies  [1978],  or  LePage  [1961].  Since  r^  and  r2  are 
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positive  integers  in  the  Erlang  distribution,  we  have  and  X2 

as  (possibly  repeated)  linear  factors  of  the  denominator. 

r,  ,  r. 


r  1 

[  ^  =  V 

ri-i 

[  s+X^  J 

[  S+X2  J  1 

S  fs+X,)^ 

4 

■2  B 


(5.29) 


To  solve  for  the  constants  A^,  we  first  rewrite  Eq  (5.29) 


21S 


-  I  S-i 


r,-i 


1=1 


(®«i)  "  1. 


2  B 

r,  r-.  r2"^ 


(5.30) 


We  can  inmediately  solve  for  Aq  by  setting  s=-X^  in  Eq  (5.30). 
Hence, 

1 


= 


(5.31) 


(^2"^l) 


Then  by  taking  successive  derivatives  of  Eq  (5.30)  with  respect 
to  s  and  evaluating  both  sides  at  s=-Xj^,  we  have 


(j) 


A.  = 
D 


j! 


ds 


TdT 


Js=-X, 


(5.32) 
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(-1)^ 


-it 


=  (-1) 


rj+j-l 


(^2-H) 


(5.33) 


for  j=l, . . .  ,rj^-l.  Note  that  in  solving  for  each  of  the  A^,  all 
terms  on  the  right  side  of  Eq  (5.30)  vanish  except  for  one. 

We  can  take  similar  steps  to  solve  for  the  B^.  We  get 


®0  = 


(j) 


B.  = 

^  j! 


ds 


w 


[  I 

L  Js=-X. 


(5.34) 


(5.35) 


(-1)’ 


j! 


^1+3 


=  (-1) 


r^+j-1 


V3 


(5.36) 


for  j=l, . . . ,r2~l. 

Hence,  we  can  write  the  Laplace  transform  of  the  p.d.f .  of 


y  as 
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(5.39) 

We  now  use  the  linearity  property  of  the  inverse  Laplace 
transform  to  invert  this  transform  term  by  term.  We  obtain  the 
p.d.f .  of  Y,  fy(y)  . 
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y  >  0 
(5.41) 

The  result  is  simply  a  sun  of  exponential  and  gamra-type  terms. 
There  are  defining  the  p.d.f.  of 

Y=X^+X^. 

The  practical  method  described  above  can  be  extended  to 

find  the  exact  distribution  of  the  svm  of  an  arbitrary  nunber 

of  independent  Erlang  randcm  variables.  Suppose  we  want  the 

n 

distribution  of  Y=  V  X.  where  X,  ,...,X  are  independent  Erlang 

i=l  ^  1 

distributed  random  variables  with  parameters  r^  and 
i=l,...,n.  Also  assume  X^^Xj  for  i#j.  If  any  scale 
parameters,  X,  were  common  among  the  random  variables  in  the 
sun,  they  could  be  immediately  combined  into  another  Erlang 
random  variable  using  the  reproductive  property. 

He  again  decompose  the  product  of  Laplace  transforms  of 
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the  individtial  p.d.f.'s  by  partial  fractions. 


(5.42) 


n  A 


y 

k 


(5.43) 


n  r  r 


(s+Xj^ 


(5.44) 


Vl  A. 

Y".  K,rj^-Tn-1 


(s+XjJ 


(5.45) 


The  constants  A.  .  can  be  found  as 

^  /  J 


(5.46) 
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for  j=0,l/ . . .  ,rj^-l  and  k=l,2,...,n. 

We  again  use  the  linearity  property  of  the  inverse  Laplace 
transform  to  invert  this  transform  term  by  term.  We  obtain  the 

p.d.f.  of  fy(y). 

V  A,  -C 

/  K,r.-m-l  m! 
mpO 

y  >  0 
(5.47) 

where  ^  are  as  given  in  Eq  (5.46)  for  j=0,l, . . .  ,rj^-l  and 
k=l,2, . . . ,n. 

Eq  (5.47)  gives  the  exact  p.d.f.  for  the  sun  of  an 

arbitrary  nunnber  of  independent  Erlamg  randan  vari^J:>les  with 

n 

different  scale  parameters  X . .  Note  that  there  are  V  r,  terms 

k=l  ^ 

in  the  sun  and  each  term  has  either  an  exponential  or  gamma- 
type  form. 

This  new  result  was  used  to  find  the  exact  distribution  of 
the  sun  of  two  independent  Erlang  variates.  Tliis  exact 
distribution  is  compared  to  the  approximate  H-function 
distribution  in  Section  5.3.2  below. 
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5.2.  MEASURES  OF  MERIT 

There  is  no  universally  accepted  meastire  to  determine  how 
"close"  an  approximation  is  to  the  exact  distribution.  It  is 
likely  that  any  standard  statistical  test  for  goodness  of  fit 
(e.g.  Chi-Square  or  Kolonogorov-Smimov)  would  fail  to  reject 
the  null  hypothesis  that  the  approximate  H-f  unction 
distribution  and  the  exact  distribution  were  equal.  However, 
this  is  because  these  tests  are  not  very  powerful  or 
discriminatory,  not  necessarily  because  the  distributions  are 
nearly  identical . 

One  way  to  cann{>are  the  approximate  H-function  distribution 
to  the  exact  distribution  involves  cori>aring  the  corresponding 
moments.  Assuning  ZSPGW  successfully  "solved"  the  system  of 
nonlinear  equations,  there  should  at  least  be  a  perfect  match 
of  the  moments  used.  For  sums  of  variates  where  an  euaalytic 
solution  for  the  H-function  parameters  was  possible,  at  least 
three  moments  should  match  perfectly.  Higher  order  moments  of 
the  approximate  H-function  distribution  may  be  in  error  to  seme 
degree. 

Although  there  is  no  standard  measure  of  the  "closeness" 
of  two  distributions,  there  are  several  measures  that  are 
comnonly  used.  Hiese  are  described  below  and  conputed  for  the 
examples  listed.  The  measures  were  included  in  a  FORTRAN 
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canriputer  program  which  compares  the  approximate  H-function 
distribution  to  the  exact  distribution  of  the  s\an. 

These  measures,  when  considered  collectively,  give  some 
idea  of  how  well  the  approximate  H-function  distribution 
matches  the  exact  distribution  of  the  sum  of  independent  H- 
f unction  variates.  Of  course,  the  mecisures  can  only  be 
computed  when  it  is  possible  to  find  the  exact  distribution  of 
the  sxxn.  In  many  cases,  the  exact  distribution  of  the  sum  of 
independent  H-fimction  variates  is  very  difficult  to  obtain. 

In  the  descriptions  that  follow,  let  the  density  of  the 
exact  distribution  of  the  sxxn  of  independent  variates  be 
represented  by  f(y),  let  the  density  of  the  approximate  H- 
function  distribution  be  H(y),  and  let  each  density  be 
evaluated  at  n  equally  spaced  values  of  the  varied^le.  The 
conrnon  interval  between  consecutive  y  values  is  Ay.  The 
densities  should  be  evaluated  over  the  whole  raoige  of  values 
the  variable  is  likely  to  assvme. 

5.2.1.  ESTIMATED  SUMS  OF  SOOftRES  OF  ERROR 

The  estimated  suns  of  squares  of  error  (SSE)  is  obtained 
by  adding  the  squared  difference  between  f(y)  and  H(y)  for  all 
n  values  of  the  variable  where  the  functions  were  evaluated. 
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Estimated  SSE  = 


J,  [  ‘W  -  “W  ]' 


(5.48) 


5.2.2.  ESTIMATED  HERN  SQUARED  ERRC» 

Unless  f(y)  and  H(y)  are  identical,  the  estimated  SSE  vd.ll 
increase  if  n  is  increased.  This  limitation  of  the  estimated 
SSE  is  corrected  by  dividing  by  n,  yielding  the  estimated  mean 
squared  error  (MSE) . 

o  1.  ■  4.  j  Estimated  SSE 

Estimated  MSE  =  - -  (5.49) 


5.2.3.  MAXIMUM  ABSOLUTE  DIFFEREMCE 

The  maximun  absolute  difference  (MAD)  measures  the  maximum 
vertical  distance  of  the  two  densities  over  the  evaluated  y 
values. 


MAD  =  max 

i=l, . . . ,n 


(5.50) 


5.2.4.  INTEGRATED  ABSOLUTE  DENSITY  DIFFERENCE 
The  integrated  absolute  density  difference  (lADD)  is  a 
measure  of  the  positive  area  between  the  two  densities.  If 
f(y)=0  for  all  ysO, 

r®  I  I 

lADD  =  I  f(y)  -  H(y)  dy  (5.51) 

Since  the  area  under  a  valid  p.d.f.  is  unity,  0  <  lADD  <  2. 
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lADD  can  be  estinated  by  siatming  the  positive  area  of 
rectangles  (or  trapezoids)  with  width  Ay  and  height 
j  f  “  ^(^ij  I  ■  cornputer  program  uses  the  larger  of  the 

vertical  distances  of  the  densities  at  the  left  and  right 
endpoints  of  each  interval  as  the  height  of  each  rectangle. 

Estimated  lADD  =  J]  ^  Ay  max(A,B)  j  (5.52) 

where 


A  = 

1  1  (5.53) 

B  = 

1 

>1 

"(='J  1 

(5.54) 

l-h 

1 

X 

1  =  0. 

Therefore,  the  estimated  lADD  will 

be  an  upper  bound  for  the  lADD. 


*5.3.  DEMC^STRATED  RESULTS 

As  shown  below,  the  new  technique  was  successful  in 
finding  an  H-function  distribution  which  closely  approximates 
the  exact  distribution  of  the  sum  of  independent  H-function 
variates . 


5.3.1.  EXAMPLE  1  -  SUM  OF  IHREE  INDEPENDEMT.  IDEWTICaiLLY 

DISTRIBUTED  Gft»1A  VARIATES 

In  this  example,  all  moments  of  the  approximate  H-function 
distribution  were  identical  to  those  of  the  exact  distribution 
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of  the  sum.  This  occurred  whenever  the  r^roductive  property 
w^ls  applicable.  In  these  cases,  the  new  technique  always  fotnd 
the  H-function  distribution  whicdi  represented  the  exact 
distribution.  The  measures  of  merit  were 
Estimated  SSE  =  0.0000000004 
Estimated  MSE  =  0.0000000000 
MAD  =  0.0000133056 

Estimated  lADD  =  0.0000050304 

As  expected,  all  measxires  of  merit  show  a  very  close  fit.  Even 
the  small  amount  of  error  is  probably  due  to  ccnputer  roimdoff 
error  in  the  evaluation  of  the  H-function  distribution  by  the 
sum  of  residues.  If  the  approximate  H-function  distribution 
were  graphically  compared  to  the  exact  distribution,  the  two 
graphs  would  be  indistinguishable. 

5.3.2.  EXAMPLE  2  -  SOM  OF  TWO  INDEPEMPBiT  ERLANG  VARIATES 

WITH  DIFFERENT  X 

Since  the  parameter  X  was  different  for  the  Erlang 
variates  in  this  example,  the  reproductive  property  of  the 
gamnna  or  Erlang  distribution  did  not  apply.  Only  the  zeroth 
throiagh  third  moments  were  identical  between  the  approximate  H- 
f unction  distribution  and  the  exact  distribution.  The  measures 
of  merit  in  this  example  were 

Estimated  SSE  =  0.0006853443 
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Estinated  MSE  =  0.0000062304 
MAD  =  0.0100130852 
Estimated  lADD  =  0.0065133019 

All  measures  of  merit  show  a  close  fit  of  the  approximate  H- 
f unction  distribution  to  the  exact  distribution.  The 
approximate  H-fmction  distribution  deviated  the  most  from  the 
exact  distribution  near  y=0.0,  achieving  its  maxinun  absolute 
deviation  of  0.01  at  y=0.1.  In  the  more  critical  right  tail, 
all  absolute  residual  valvies  were  less  than  0.0006. 

The  two  distributions  are  compared  graphically  in  Figure  3 
below.  The  densities  are  nearly  identical. 
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continuous  at  y=l.  In  general,  the  new  H-f unction  technique 
was  less  effective  when  the  derivative  of  the  exact 
distribution  was  not  contintious. 

Of  course,  either  linear  function  in  Eq  (5.13)  could  be 
exactly  represented  as  an  H-fmcticn.  With  a  shift  in  the 
argunnent  of  one  H-function,  it  is  even  possible  to  perfectly 
represent  the  exact  distribution  of  the  stm  as  a  sxm  of  H- 
f unctions.  But  in  the  new  tedtoique,  we  try  to  fit  a  single  H- 
function  to  the  moments  of  the  svm.  This  produces  sm  H- 
function  variate  which  (correctly)  cannot  assunne  values  less 
than  zero  or  greater  than  two.  Even  the  moments  of  the 
approximate  H-function  distribution  were  reasonably  close  to 
those  of  the  exact  distribution.  Still,  the  measures  of  merit 
were 

Estimated  SSE  =  0.0920133337 
Estimated  MSE  =  0.0023003333 
MAD  =  0.1511936368 
Estimated  lADD  =  0.0752034552 

The  approximate  H-function  distribution  was  not  as  effective  in 
representing  the  exact  distribution  as  in  the  earlier  exanples. 
The  approximate  H-function  distribution  achieved  its  maximum 
absolute  deviation  of  0.15  at  y=1.0.  It  more  closely 
approximated  the  exact  distribution  in  both  tails.  The  two 


distributions  are  contpared  graphically  in  Figure  4  below. 
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fit  is  obtained,  at  least  over  the  range  of  the  variable  from 
zero  to  two.  The  meausures  of  merit  were 
Estimated  SSE  =  0.0620352924 
Estimated  MSE  =  0.0015508823 
MAD  =  0.1107462928 

Estimated  lADD  =  0.0581953790 

The  drawback  of  this  slightly  better  fit  is  that  this  H- 
ftmction  distribution  is  positive  for  all  y>0.  At  y=2.0,  the 
first-order  H-function  distribution  has  f(2)=.0348  and 
F(2)=.9966.  Therefore,  the  probability  of  achieving  a  value 
greater  than  two  with  the  first-order  H-function  distribution 
is  .0034.  Of  course,  this  event  is  impossible  and  the 
probability  should  be  zero.  Most  of  the  measures  given  above 
would  increase  (worsen)  if  they  were  computed  over  a  larger 
range  of  the  variable  since  the  exact  distribution  has  f(y)=0 
for  y>2  but  the  first-order  H-function  distribution  has  f(y)>0 
for  y>2. 

5.3.4.  EXAMPLE  4  -  SOM  OF  TWO  INraTEMPaiT.  IDEWTICALLY 

DISroiBUTB)  BETA  VARIATES 

This  exanple  demonstrates  the  technique  for  the  svm  of  two 
Type  VI  convergent  H^  ^  variates  over  (0,1)  where  the  first 
derivative  of  the  exact  distribution  is  continuous  at  y=l.  The 
exact  distribution  is  given  by  Eq  (5.17).  Again,  the  moments 
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of  the  approximate  H-f motion  distribution  were  reasonably 
close  to  those  of  the  exact  distribution.  The  mecisures  of 
merit  were 

Estimated  SSE  =  0.0358601689 
Estimated  MSE  =  0.0008965042 
MM)  =  0.0579816425 
Estimated  IM)D  =  0.0442126253 

The  approximate  H-f motion  distribution  more  closely  matched 
the  exact  distribution  than  in  the  previom  example,  but  not  as 
close  as  in  the  gamma  and  Erlang  examples.  While  the  two 
fmctions  which  represent  the  exact  distribution  have  a 
continuous  first  derivative  at  y=l,  the  second  derivative  is 
not  continuous  there.  The  approximate  H-f motion  distribution 
has  continuous  derivatives  of  all  orders  for  y  €  (0,2).  The 
two  distributions  are  compared  graphically  in  Figure  5  below. 


CHAPTER  6 


CONCLUSIONS  AND  REOOMEKDATIONS  FOR  FURTHER  STUDY 


This  thesis  has  presented  a  practical  method  to  find  an  H- 
function  distribution  for  the  sun  of  two  or  more  independent  H- 
f unction  variates.  The  method  finds  the  mcments  about  the 
origin  of  each  variate  in  the  s\jrn  and  uses  these  to  find  the 
mcments  of  the  sun.  The  parameters  of  an  H-function  are  then 
estimated  using  the  method  of  moments.  This  produces  an  H- 
f unction  distribution  whose  moments  closely  approximate  those 
of  the  sum  of  independent  H-fvmcticjn  variates. 

The  new  technique  is  especially  useful  because  the  H- 
f unction  c:an  exacrtly  represent  nearly  every  named  statistical 
distribution  of  positive  random  variables.  By  rewriting  the 
distributions  of  the  random  variables  in  the  sun  in  their  H- 
function  representations,  an  analyst  can  find  the  sun  of  two  or 
more  independent  variates  with  practically  any  distributional 
form. 

Further,  simple  formulas  exist  which  immediately  give  the 
distribution,  as  an  H-function,  of  the  random  variable  defined 
as  the  prcxiuct,  quotient,  or  power  of  independent  H-function 
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variates.  The  new  technique  for  the  sum  allows  an  analyst  to 
exploit  these  powerful  H-ftmction  results  to  find  the  H- 
f unction  distribution  of  more  complicated  algebraic 
ccnbinations  of  independent  random  variables. 

When  the  exact  distribution  of  an  algebraic  combination  of 
independent  random  variables  is  difficult  or  iitpossible  to 
obtain,  an  analyst  may  resort  to  computer  simulation  to  analyze 
its  distributional  properties.  The  method  presented  here,  when 
combined  with  the  other  H-function  results,  give  an  H-ftmction 
distribution  for  the  algebraic  ccmbination,  possibly  precluding 
the  need  to  rely  on  computer  simulation. 

A  FORTRAN  computer  program  which  inplanents  the  new 
technique  is  demonstrated  throxigh  five  exanples.  Of  course, 
the  program  was  also  tested  on  many  other  exanples.  When  it 
was  possible  to  find  the  exact  distribution  for  the  sum  of 
independent  random  variables,  the  approximate  H-function 
distribution  was  compared  to  the  exact  distribution. 

In  cases  vAiere  the  exact  distribution  of  the  svm  is  widely 
known  (e.g.  reproductive  distributions),  the  method  always 
found  the  correct  H-function  representation  for  the  sum.  In 
other  cases,  only  the  moments  used  to  estimate  the  H-function 
parameters  agreed  with  those  of  the  sum.  Other,  higher-order 
moments  were  in  error  to  some  degree.  The  method  was  less 
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effective  in  finding  an  H-fvmction  distribution  which  closely 
approxinated  the  exact  distribution  of  the  stm  of  independent 
variates  over  a  restricted  range.  Itiese  sums  produce  an  exact 
distribution  with  two  functional  fonns  over  distinct  ranges  of 
the  variable  and  do  not  have  continuous  derivatives  of  all 
orders  at  y=l. 

The  exact  distribution  of  the  general  sum  of  independent 
Erlang  variates  with  different  scale  parameters,  X,  was  derived 
in  Section  5.1.3.  An  Erlang  variate  is  sinply  a  gamma  variate 
with  an  integer  shape  parameter  r.  The  derivation  used  partial 
fractions  to  decompose  the  product  of  Laplace  transforms  of  the 
individual  densities.  This  produced  a  sm  of  terms,  each  of 
which  could  easily  be  inverted  from  transform  space,  yielding 
the  desired  density  of  the  sum  of  independent  variates. 

Throughout  the  thesis,  a  nunber  of  other  new  results 
relating  to  the  H-function  were  given.  A  previously  unstated 
restriction  on  the  variable  in  the  H-function  representations 
of  power  functions  and  beta-type  fimctions  was  highlighted. 
Several  ways  of  overcoming  this  limitation  when  representing 
mathematical  functions  were  presented.  These  include  scaling 
the  variable,  allowing  a  slightly  different  function  to  be 
represented  for  a  range  of  the  variable,  or  the  use  of  complex 


parameters . 
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The  restric±ion,  however,  is  an  advantage  when 
representing  certain  statistical  distributions.  For  randcm 
variables  with  a  restricted  range  of  the  variable,  the  H- 
f unction  representation  eliminates  the  need  to  specify  the 
range  for  which  the  density  is  nonzero.  The  H-f unction 
dis;.  'ibution  exactly  represents  the  desired  density  over  the 
appropriate  range  and  is  zero  elsewhere. 

Through  scaling  the  variable  of  these  H-f unctions  with  the 
parameter  c,  the  value  where  the  H-functim  changes  from 
representing  the  special  case  to  taking  the  value  zero  can  be 
changed.  This  allowed  the  H-fvjnction  to  represent  the  power 
fxmctions  and  beta-type  functions  over  a  more  general  range  arid 
also  gave  the  H-function  representation  of  the  Unit  Step 
Function  and  its  conpl ament.  The  same  technique  was  used  to 
give  the  H-function  representation  of  the  Three-Parameter  Beta 
p.d.f.,  the  Power  Function  p.d.f.  over  (0,M),  the  Uniform 
p.d.f.  over  (0,M),  and  the  Pareto  p.d.f.  over  (e,oo). 

Analytic  solutions  were  derived  for  the  system  of 
nonlinear  equations  generated  by  the  method  of  moments  to 
estimate  the  H-function  parameters  for  certain  distributions. 
These  distributions  have  restricted  ranges  of  the  variable  and 
must  be  represented  by  Type  VI  convergent  H-f unctions. 

Many  new  H-function  representations  of  other  mathematical 
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functions  were  also  given,  including  those  of  the  carpi ementary 
error  function,  erfc  x,  the  caiplementary  incorplete  ganma 
function,  r(a,x),  the  caiplementary  incorplete  beta  function, 
B^(a,^),  and  the  inverse  hyperbolic  cotangent  function, 
arccoth  z. 

Four  new  theorems  show  when  and  how  a  generalizing 
constant  may  be  present  in  an  H-f unction  representation.  Many 
generalized  H-f unction  representations  are  given,  including 
those  of  every  cunulative  distribution  function  of  an  H- 
function  variate.  The  generalizing  constant  was  also  possible 
in  the  H-f unction  representations  of  power  functions,  the  error 
f motion  and  its  carpi ement,  the  incorplete  gamma  function  and 
its  carpi  ement,  the  incorplete  beta  function  and  its 
carpi  ement,  many  inverse  trigonanetric  and  hyperbolic 
functions,  and  the  logarithmic  functions. 

Sufficient  convergence  COTiditions  were  developed  for  the 
alternate  definition  of  the  H-function.  These  show  where  the 
H-function  may  be  evaluated  by  the  sun  of  LHP  residues  or  the 
negative  sun  of  RHP  residues,  without  first  changing  the  form 
of  the  alternate  definition  of  the  H-function  to  that  of  the 
primary  definition. 

The  hierarchical  strxicture  among  classes  of  H-functions 
was  given  through  seven  new  theorems.  Every  cl2L5s  of  H- 
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functions  is  wholly  contained  in  many  higher-order  classes  of 
H-functions  through  the  application  of  the  dv^jlication, 
triplication,  and  multiplication  formulas  for  the  ganim 
function.  In  the  other  direction,  the  hierarchical  structure 
gives  new  reduction  properties  that  are  less  obvious  than  the 
known  reduction  properties. 

Several  new  results  were  obtained  for  powers  of  first- 

order  H-f unction  variates.  The  square  of  a  random  variable 

with  a  Rayleigh  distribution  with  parameter  X  has  an 

exponential  distribution  with  parameter  X.  If  X  has  a  Weibull 

distribution  with  parameters  ^  and  X,  then  y=X^  where  P>0  has  a 

Weibull  distribution  with  parameters  and  X.  The  square  of 

a  random  variable  with  a  half-normal  distribution  with 

parameter  a  has  a  ganrna  distribution  with  parameters  r  = 

and  X  =  — —  The  square  of  a  random  variable  with  a  Maxwell 
2a^ 

distribution  with  parameter  6  has  a  gamma  distribution  with 
3  1 

parameters  r  =  — ^  and  X  =  — « —  The  square  of  a  random 

variable  with  a  Maxwell  distribution  with  parameter  6-/7  heis 
a  Chi-Square  distribution  with  parameter  p=3.  The  square  of  a 
random  variable  with  a  Rayleigh  distribution  with  parameter 
X  =  — ^  has  a  Chi-Square  distribution  with  parameter  t^=2. 

Many  areas  of  further  research  proposed  by  Eldred  [1979, 
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pp.  258-9],  Cook  [1981,  pp.  145-147],  Kellogg  [1984,  pp.l70- 
171],  and  Jacobs  [1986,  pp.  149-152]  remain  unfinished. 
Readers  interested  in  H-function  research  should  consider  those 
reccnmendatians  and  some  additional  areas  that  follow. 

In  Bayesian  statistics,  the  parameters  of  a  distribution 
are  considered  to  be  random  variables  with  a  particular 
distribution.  After  sampling  from  the  distribution,  the 
information  gained  about  the  parameters  is  reflected  by 
updating  the  prior  distribution  and  creating  a  posterior 
distribution.  To  my  knowledge,  a  Bayesian  approach  has  not 
been  attempted  with  the  H-function. 

Nxmerical  methods  exist  for  certain  cormonly-used 
transfonns  such  as  the  Laplace  and  Fourier  transforms.  If 
these  methods  are  sufficiently  accurate  and  will  work  in  both 
directions  (finding  a  transform  and  inverting  a  transform),  it 
would  be  possible  to  find  the  distribution  (in  tabuleu:  form)  of 
sixns  or  differences  of  independent  random  variables.  If  a 
nvmerical  method  to  invert  a  Mellin  transform  exists,  it  could 
be  used  to  evaluate  the  H-function. 

The  computer  program  by  Cook  [1981,  Cook  and  Barnes,  1981] 
will  evaluate  an  H-function  by  the  sum  of  residues  for  certain 
values  of  the  real  variable,  x,  when  given  the  parameters  of 
the  H-function.  The  program's  capabilities  should  be  extended 
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to  allow  corplex  parameters  of  the  H-fmction  and  to  evaluate 
the  H- function  at  values  of  the  ccnnplex  variable,  z. 

The  nixnerical  solution  of  the  system  of  nonlinear 
equations  with  ZSPOH  is  scmetimes  nixnerically  vnstable. 
Several  other  nunerical  methods  were  tried,  but  none  worked  as 
well  eis  ZSPON  in  solving  the  type  of  nonlinear  equations 
generated  \ising  the  method  of  moments  with  the  H-f unction. 
Further  research  might  develop  a  better  way  to  solve  the  system 
of  nonlinear  equations. 

While  the  H-f  unction  has  been  shown  to  exactly  represent 
many  mathematical  functions  and  statistical  distributions, 
there  are  a  few  functions  and  distributions  vdiich  have  not  yet 
been  shown  to  be  H-f  unctions.  Conversely,  no  one  has  proven 
that  these  functions  cannot  be  represented  as  H-f  motions. 
Exairples  include  the  ganma  function,  certain  trigonometric  and 
hyperbolic  functions,  and  the  Lognormal  and  Logistic  densities. 
A  related  issue  is  to  prove  or  disprove  that  the  product  of  two 
H-f  unctions  is  another  H-f  unction. 

To  my  knowledge,  the  H-f  unction  has  not  been  applied  to 
the  study  of  complex-valued  random  variables.  The  H-fmction 
seems  the  appropriate  tool  to  analyze  random  variables  which 
take  values  in  the  complex  plane. 

The  relationship  between  the  H-f  unction  and  the  infinite 
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sunnrability  of  integrals  described  in  Section  2.7  needs  to  be 
more  extensively  stixlied.  Is  there  a  particular  sunnrability 
scheme  the  H-f unction  always  follows  to  obtain  the  moments  of 
the  functions  representable  as  H-f unctions? 

Theorems  similar  to  those  in  Section  2.4.6  could  be 
developed  to  show  \;TOrade  and  reduction  results  for  H-f  unctions 
by  applying  the  triplication  and  multiplication  formulas  for 
the  gainra  function.  Then,  the  search  could  begin  for  a  class 
of  H-f unctions  which  includes  most  of  the  known  special  cases. 
For  example.  Figure  1  showed  the  H^  ^  class  includes  many 
common  statistical  distributions,  including  all  those  in  the 
H^  J  class.  The  half-Cauchy,  half -Student,  and  F 
distributions,  however,  are  in  the  H^  J  class  of  H-f  unction 
distributions.  Using  combinations  of  theorems  based  on  the 
duplication,  triplication,  or  multiplication  formulas,  is  there 
a  class  of  H-functions  which  includes  both  the  H^  ^  and  H^  ^ 
classes  of  H-functions? 

The  H-f unction  continues  to  be  an  extremely  rewarding  area 
of  research  for  me.  I  am  convinced  that  we  have  much  to 
discover  about  its  ability  to  simplify  many  difficult  problems 
in  mathematics,  probability,  and  statistics.  Those  who  pursue 
H-fxmction  research  should  find  it  both  challenging  and 
rewarding . 


APPEJOIX  A 


OUTPUT  ETROM  COMPUTER  PROGRftMS 


Listed  below  are  the  ccnplete  output  files  from  the 
FCXnRAN  corputer  programs  for  each  of  the  five  examples  in 
Section  4.4.  The  output  from  the  program  which  finds  a  H- 
function  distribution  for  the  sun  of  independent  H-function 
variates  is  followed  by  the  output  from  a  different  program 
[Cook,  1981;  Cook  and  Barnes,  1981]  vMch  evaluates  the  p.d.f. 
and  c.d.f.  of  the  resulting  H~function  distribution.  For  the 
first  four  examples  for  which  it  was  possible  to  find  the  exact 
distribution,  the  output  from  the  program  which  compares  the 
approximate  H-fxmction  distribution  to  the  exact  distribution 
is  also  presented. 

EataMPLE  1  ^  ^  OT  THREE  INTgRERIgRT.  IDERTICALLY 
DISTRIBUTED  GR>»iA  VARIATES 


PROGRAM  SUMVAR  RUN  IN  DOUBLE  PRECISION 
FOR  THE  SUM  OF  3  INDEPENDENT  H-EDNCmON  VARIATES 
INPUT  PARAMETERS  FCXl  VARIABLE  NUMBER  1  THE  SUM: 
M,  N,  P,  Q  =  10  0  1 

K,  C  =  2.0000000  2.0000000 

b,  B  =  1.00000  1.00000 

THE  MOMENTS  ABOUT  THE  CXIIGIN  ARE: 

THE  ZEROTH  MOMENT  IS  1.0000000 

THE  FIRST  MOMENT  IS  1.0000000 
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THE  SEOaiD  MCMEDT  IS  1.5000000 

IHE  THIRD  MQMEirr  IS  3.0000000 

TOE  POURTO  MCMEJIT  IS  7.5000000 

TOE  FIFTO  MCMEWT  IS  22.5000000 

INPUT  PARAMETERS  FOR  VARIABLE  NUMBER  2  OF  TOE  SUM 
M,  N,  P,  Q  =  1001 

K,  C  =  2.0000000  2.0000000 

b,  B  =  1.00000  1.00000 

TOE  MCMERTS  ABOUT  TOE  ORIGIN  ARE: 

TOE  ZroOTO  MOMEWT  IS  1.0000000 

TOE  FIRST  MCMERT  IS  1.0000000 

TOE  SBOOID  MOMENT  IS  1.5000000 

THE  TOIRD  MCMERT  IS  3.0000000 

TOE  FOURTH  MQMERT  IS  7.5000000 

THE  FIFTO  MQMERT  IS  22.5000000 

INPUT  PARAMETTERS  FC»  VARIABLE  NUMBER  3  OP  TOE  SUM 
M,  N,  P,  Q  =  1001 

K,  C  =  2.0000000  2.0000000 

b,  B  =  1.00000  1.00000 

TOE  MDMERTS  ABOUT  THE  ORIGIN  ARE: 

THE  ZEROTH  MCMERT  IS  1.0000000 

THE  FIRST  MCMERT  IS  1.0000000 

THE  SBOCXO}  MCMERT  IS  1.5000000 

THE  TOIRD  MCMERT  IS  3.0000000 

THE  FOURTH  MCMERT  IS  7.5000000 

THE  FIFTO  MCMERT  IS  22.5000000 

FOR  THE  RANDOM  VARIABLE  GIVER  AS  THE  SUM 
OF  THE  ABOVE  INDEPERE®RT  H-FUNCmON 
VARIATES,  THE  MCMERTS  ABOUT  THE  ORIGIN  ARE: 

THE  ZEROTH  MCMERT  IS  1.0000000 

THE  FIRST  MOMENT  IS  3.0000000 

THE  SECOND  MOMENT  IS  10.5000000 

THE  TOIRD  MQMERT  IS  42.0000000 

THE  FOURTH  MCMERT  IS  189.0000000 

THE  FIFTH  MCMERT  IS  945.0000000 

PROGRAM  H_FIT  RUN  WITH  DOUBLE  PRECISICW 

INPUT  WAS  OF  TYPED 

DEFAULT  INITIAL  GUESS  WAS  USED 

THE  ZBK3TH  MOMENT  WAS  USED  IN  THE  FIT 

RESULTS  OF  ZSPOW  - 

NUMERATOR: 

SMALLB(1)=  4.999999888533181114 


BIGB(1)=  0.999999992560842613 

imi4IIATQR: 

VMiUES  OF  K  &  C  ARE: 

K=  0 . 016666666666666667 

C=  2 . 000000000000000000 

F1K]RM=  0.000000000000000000 

ESTIMATED  H-FUNCrTIGN  PARAMETE3{S 
FCR  THE  SUM  OF  THE  INDEPEafDEHT 
H-FHNCTION  VARIATES  ARE: 

M,  N,  P,  Q  =  10  0  1 

K,  C  =  0.0166667  2.0000000 

b,  B  =  5.00000  1.00000 

THE  MCMEHTS  ABOUT  THE  C»IGIN  OF  THIS 
H-ETJNC?riON  ARE: 

THE  ZEROTH  MOMENT  IS  1.0000000 

THE  FIRST  MOMENT  IS  3.0000000 

THE  SBOCXID  t^OMENT  IS  10.5000000 

THE  THIRD  MOMENT  IS  42.0000000 

THE  FOURTH  MOMENT  IS  189.0000000 

THE  FIFTH  MOMENT  IS  945.0000000 

DETERMINE  P.D.F.(Z)  AND  C.D.F.(Z) 

FOR  VALUES  OF  Z  EROM  0.1000  TO  7.0000 

WITH  STEP  SIZE  0.1000 

FOR  THE  SIM  CF  1  TERMS,  WHERE 

THE  MAXIMUM  NUMBER  CF  POLES  TO  BE  EVALUATED  IS  100. 

CRDMP  PARAMETERS:  NUMBER  OF  (XMPLEX  VALUES  =  1001, 

PERCENT  OF  HIGHEST  Z  VALUE  =  1.00,  AXIS  POINT  A  =  1.3653 

FORM  PC»  OVERALL  PROBLEM  (WHERE  YJ  =  XJ**PJ) : 

Z  =  Y1 

VARIATE  X  1  IS  TYPE  NUMBER  4 

INPUT  PARAMEmiS  ARE  THETA  =  0.01667,  PHI  =  2.00000 

AND  POWER  =  1.00000 

THE  P.D.F.  FOR  VARIATE  X  1  IS  GIVEN  BY: 

1  0 

0.01667  H  (  2.00000  X): 

0  1  (  5.000,  1.000) 

THE  P.D.F.  FOR  TERM  1  OF  THE  SUM  IS  GIVB<  BY: 

1  0 

0.01667  H  (  2.00000  Z),  WHERE 

0  1 

(BA(I),GBA(I));  (  5.000,  1.000)  ( 

CONVERGENCE  TYPE  =  1 

D  =  1.00  E  =  -1.00  L  = 


4.50 


R  =  1.0000 
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4.6000 

0.110989 

0.895926 

4.7000 

0.101186 

0.906529 

4.8000 

0.092040 

0.916185 

4.9000 

0.083540 

0.924959 

5.0000 

0.075667 

0.932914 

5.1000 

0.068398 

0.940112 

5.2000 

0.061710 

0.946613 

5.3000 

0.055572 

0.952472 

5.4000 

0.049956 

0.957745 

5.5000 

0.044830 

0.962480 

5.6000 

0.040164 

0.966726 

5.7000 

0.035926 

0.970527 

5.8000 

0.032086 

0.973925 

5.9000 

0.028614 

0.976957 

6.0000 

0.025481 

0.979659 

6.1000 

0.022659 

0.982063 

6.2000 

0.020123 

0.984200 

6.3000 

0.017847 

0.986096 

6.4000 

0.015808 

0.987777 

6.5000 

0.013985 

0.989265 

6.6000 

0.012357 

0.990581 

6.7000 

0.010905 

0.991742 

6.8000 

0.009612 

0.992767 

6.9000 

0.008462 

0.993669 

7.0000 

0.007440 

0.994464 

NUMBER  OF  POLES  EVALUATED  =  83 


FBOGRAM  OCMPAR  RUN  IN  DOUBLE  PRECISION 


X 

Y  (H-FCN) 

Y  (EVAL) 

RESIDUAL 

0.1000 

0.0000043666 

0.0000043666 

0.0000000000 

0.2000 

0.0001144013 

0.0001144013 

0.0000000000 

0.3000 

0.0007112599 

0.0007112599 

0.0000000000 

0.4000 

0.0024539353 

0.0024539352 

0.0000000001 

0.5000 

0.0061313241 

0.0061313240 

0.0000000001 

0.6000 

0.0124911263 

0.0124911264 

-0.0000000001 

0.7000 

0.0221042938 

0.0221042943 

-0.0000000005 

0.8000 

0.0352839727 

0.0352839739 

-0.0000000012 

0.9000 

0.0520572460 

0.0520572483 

-0.0000000023 

1.0000 

0.0721788138 

0.0721788177 

-0.0000000039 

1.1000 

0.0951731110 

0.0951731171 

-0.0000000061 

1.2000 

0.1203921499 

0.1203921587 

-0.0000000088 

1.3000 

0.1470787062 

0.1470787183 

-0.0000000121 

1.4000 

0.1744272435 

0.1744272593 

-0.0000000158 

1.5000 

0.2016376070 

0.2016376269 

-0.0000000199 

1.6000 

0.2279587426 

0.2279587669 

-0.0000000243 

1.7000 

0.2527214162 

0.2527214452 

-0.0000000290 

198 


1.8000  0 

1.9000  0 

2.0000  0 

2.1000  0 

2.2000  0 

2.3000  0 

2.4000  0 

2.5000  0 

2.6000  0 

2.7000  0 

2.8000  0 

2.9000  0 

3.0000  0 

3.1000  0 

3.2000  0 

3.3000  0 

3.4000  0 

3.5000  0 

3.6000  0 

3.7000  0 

3.8000  0 

3.9000  0 

4.0000  0 

4.1000  0 

4.2000  0 

4.3000  0 

4.4000  0 

4.5000  0 

4.6000  0 

4.7000  0 

4.8000  0 

4.9000  0 

5.0000  0 

5.1000  0 

5.2000  0 

5.3000  0 

5.4000  0 

5.5000  0 

5.6000  0 

5.7000  0 

5.8000  0 

5.9000  0 

6.0000  0 

6.1000  0 

6.2000  0 

6.3000  0 


.2753601348 

.2954252728 

.3125868609 

.3266316878 

.3374553719 

.3450509461 

.3494953096 

.3509346792 

.3495699434 

.3456426020 

.3394217767 

.3311926087 

.3212462161 

.3098712736 

.2973471912 

.2839388042 

.2698924499 

.2554332757 

.2407636152 

.2260622642 

.2114844971 

.1971626696 

.1832072739 

.1697083224 

.1567369566 

.1443471907 

.1325777187 

.1214537217 

.1109886370 

.1011858411 

.0920402391 

.0835397256 

.0756665181 

.0683983578 

.0617095716 

.0555720188 

.0499558733 

.0448303394 

.0401642223 

.0359264174 

.0320862857 

.0286140192 

.0254807806 

.0226590473 

.0201226143 

.0178466675 


0.2753601684 

0.2954253111 

0.3125869037 

0.3266317349 

0.3374554229 

0.3450510007 

0.3494953673 

0.3509347395 

0.3495700059 

0.3456426662 

0.3394218420 

0.3311926746 

0.3212462821 

0.3098713394 

0.2973472562 

0.2839388681 

0.2698925124 

0.2554333366 

0.2407636741 

0.2260623210 

0.2114845515 

0.1971627215 

0.1832073232 

0.1697083691 

0.1567370007 

0.1443472324 

0.1325777577 

0.1214537587 

0.1109886714 

0.1011858741 

0.0920402710 

0.0835397566 

0.0756665496 

0.0683983922 

0.0617096146 

0.0555720692 

0.0499559385 

0.0448304269 

0.0401643429 

0.0359265830 

0.0320865265 

0.0286143510 

0,0254812775 

0.0226597557 

0.0201235965 

0.0178480594 


-0.0000000336 

-0.0000000383 

-0.0000000428 

-0.0000000471 

-0.0000000510 

-0.0000000546 

-0.0000000577 

-0.0000000603 

-0.0000000625 

-0.0000000642 

-0.0000000653 

-0.0000000659 

-0.0000000660 

-0.0000000658 

-0.0000000650 

-0.0000000639 

-0.0000000625 

-0.0000000609 

-0.0000000589 

-0.0000000568 

-0.0000000544 

-0.0000000519 

-0.0000000493 

-0.0000000467 

-0.0000000441 

-0.0000000417 

-0.0000000390 

-0.0000000370 

-0.0000000344 

-0.0000000330 

-0.0000000319 

-0.0000000310 

-0.0000000315 

-0.0000000344 

-0.0000000430 

-0.0000000504 

-0.0000000652 

-0.0000000875 

-0.0000001206 

-0.0000001656 

-0.0000002408 

-0.0000003318 

-0.0000004969 

-0.0000007084 

-0.0000009822 

-0.0000013919 
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6.4000  0.0158078894  0.0158099039  -0.0000020145 

6.5000  0.0139846566  0.0139874081  -0.0000027515 

6.6000  0.0123566022  0.0123603635  -0.0000037613 

6.7000  0.0109048623  0.0109100486  -0.0000051863 

6.8000  0.0096119998  0.0096191863  -0.0000071865 

6.9000  0.0084623041  0.0084718905  -0.0000095864 

7.0000  0.0074402961  0.0074536017  -0.0000133056 

ESTIMATED  SSE  IS  0.0000000004 

ESTIMATED  MSE  IS  0.0000000000 

110:  MAXIMUM  ABSOLUTE  DIFFERE34CE  BE7IVIEE3<  THE 
H-ETJNCTION  AND  THE  ACTUAL  FUNCTION  IS  0.0000133056 

THE  INTEGRATED  ABSOLUTE  DENSITY  DIETTKENCE 
(AN  ESTIMATE  OF  THE  AREA  BETWEEN  THE 

H-FONCriO*  AND  THE  AdUAL  DENSITY)  IS  0.0000050304 


EXflMPLE  2  -  ^  OF  TOO  INDEPENDENT  ERLANG  VARIATES 
WITH  DIE'FERENT  X 


PROGRAM  SUMVAR  RUN  IN  DOUBLE  PRECISION 
FOR  THE  SUM  OF  2  INDEPENDENT  H-EUNCTION  VARIATES 
INPUT  PARAMETERS  FOR  VARIABLE  NUMBER  1  OF  THE  SUM: 
M,  N,  P,  Q  =  10  0  1 

K,  C  =  4.0000000  4.0000000 

b,  B  =  1.00000  1.00000 

THE  MCMENTS  ABOUT  THE  ORIGIN  ARE: 

THE  ZEROTH  MOMENT  IS  1.0000000 

THE  FIRST  MOMENT  IS  0.5000000 

THE  SECOND  MOMENT  IS  0.3750000 

THE  THIRD  MOMENT  IS  0.3750000 

THE  FOURTH  MOMENT  IS  0.4687500 

THE  FIFTH  MOMENT  IS  0.7031250 

THE  NEXT  MOMENT  IS  1.2304688 

THE  NEXT  MOMENT  IS  2.4609375 

INPUT  PARAMETERS  FCX?  VARIABLE  NUMBER  2  OF  THE  SUM: 
M,  N,  P,  Q  =  10  0  1 

K,  C  =  2.0000000  2.0000000 

b,  B  =  0.00000  1.00000 

THE  MCMENTS  ABOUT  THE  CX^IGIN  ARE: 

THE  ZEROTH  MOMENT  IS  1.0000000 

THE  FIRST  MOMENT  IS  0.5000000 

THE  SEXXyUD  MOMENT  IS  0.5000000 

THE  THIRD  MOMENT  IS  0.7500000 


THE  FOORTO  MCMEWT  IS 


TOE  FIPTO 
TOE  NEXT 
TOE  NEXT 


MOMENT  IS 
MOMENT  IS 
MOMENT  IS 


1.5000000 

3,7500000 

11.2500000 

39.3750000 


FOR  TOE  RANDOM  VARIABLE  GIVEN  AS  TOE  SOM 
OF  TOE  ABOVE  INDEPENDENT  H-FUNCTIC»I 
VARIATES,  TOE  MOMENTS  ABOUT  TOE  CMIGIN  ARE: 


TOE 

ZEROTO  MOMENT 

IS 

1.0000000 

TOE 

FIRST 

MOMENT 

IS 

1.0000000 

TOE 

SEOCM)  MOMENT 

IS 

1.3750000 

TOE 

TOIRD 

MOMENT 

IS 

2.4375000 

TOE 

FOURTH  MOMENT 

IS 

5.3437500 

TOE 

FIFTH 

MOMENT 

IS 

14.0625000 

TOE 

NEXT 

MOMENT 

IS 

43.4179688 

THE 

NEXT 

MOMENT 

IS 

154.4238281 

PROGRAM  H_FIT  RUN  WITO  DOUBLE  PRECISIC»I 

INPUT  WAS  OF  TYPEO 

DEFAULT  INITIAL  GUESS  WAS  USED 

THE  ZEROTO  MOMENT  WAS  USED  IN  TOE  FIT 

RESULTS  OF  ZSPGW  - 

NUMERATOR: 

SMALLB(1)=  3.024855875185841148 

BIGB(1)=  1.263891763119380202 

DENOMINATCR: 

VALUES  OF  K  S  C  ARE: 

K=  0 . 748958966084013356 

C=  6.533302402015585986 

FNORM=  0.000000000000000000 


ESTIMATED  H-FUNCTIOJ  PARAMETERS 
FCR  TOE  SUM  OF  TOE  INDEPENDENT 
H-ETJNCTICW  VARIATES  ARE: 


M,  N,  P,  Q  =  10 

0 

1 

K,  C  =  0.7489590 

6 

.5333024 

b,  B  =  3.02486 

1.26389 

THE  MOMENTS  ABOUT  THE  ORIGIN  OF  THIS 
H-FUNCTICR  ARE: 

THE  ZEROTO  MOMENT 

IS 

1.0000000 

THE  FIRST  MOMENT 

IS 

1.0000000 

THE  SECOUD  MOMENT 

IS 

1.3750000 

THE  TOIRD  MOMENT 

IS 

2.4375000 

THE  FOURTH  MOMENT 

IS 

5.3382458 

THE  FIFTH  MOMENT 

IS 

14.0112602 

THE  NEXT  MOMENT 

IS 

43.0835662 
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THE  NEXT  MOMEKT  IS  152.4885870 

DETERMINE  P.D.F.(Z)  AND  C.D.F.(Z) 

PC»  VALUES  Cff’  Z  FRCM  0.0500  TO  5.5000 

WITH  STEP  SIZE  0.0500 

FOR  THE  SUM  OF  1  TERMS,  WHERE 

THE  MAXIMUM  NUMBER  OF  POLES  TO  BE  EVALUATED  IS  100 . 

CRUMP  PARAMETERS:  NUMBER  OF  CXWPLEX  VALUES  =  1001, 

PERCENT  OF  HIGHEST  Z  VALUE  =  1.00,  AXIS  POINT  A  =  1.7376 

FORM  FC»  OVERALL  PROBLEM  (WHERE  YJ  =  XJ**PJ) : 

Z  =  Y1 

VARIATE  X  1  IS  TYPE  NUMBER  4 

INPUT  PARAMETERS  ARE  THETA  =  0.74896,  PHI  =  6.53330 

AND  POWER  =  1.00000 

THE  P.D.F.  FOR  VARIATE  X  1  IS  GIVM  BY: 

1  0 

0.74896  H  (  6.53330  X): 

0  1  (  3.025,  1.264) 

THE  P.D.F.  PC»  TERM  1  OF  THE  SUM  IS  GIVER  BY: 

1  0 

0.74896  H  (  6.53330  Z),  WHERE 

0  1 

(BA(I),GBA(I)):  (  3.025,  1.264)  ( 

CXRVERGERCE  TYPE  =  1 

D  =  1.26  E  =  -1.26  L  =  2.52  R  =  0.7438 

Z  PDF(Z)  CDF(Z) 

0.0500  0.026956  0.000430 

0.1000  0.104760  0.003557 

0.1500  0.211027  0.011377 

0.2000  0.326668  0.024812 

0.2500  0.438949  0.043985 

0 . 3000  0 . 540207  0 . 068521 

0.3500  0.626386  0.097754 

0.4000  0.695859  0.130880 

0.4500  0.748580  0.167060 

0.5000  0.785467  0.205474 

0.5500  0.807989  0.245367 

0.6000  0.817884  0.286063 

0.6500  0.816969  0.326975 

0.7000  0.807023  0.367609 

0.7500  0.789713  0.407555 

0 . 8000  0 . 7 66559  0 . 446483 

0.8500  0.738910  0.484136 

0.9000  0.707945  0.520319 

0.9500  0.674675  0.554892 

1.0000  0.639955  0.587762 
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1.0500 

1.1000 

1.1500 
1.2000 

1.2500 

1.3000 

1.3500 

1.4000 

1.4500 

1.5000 

1.5500 
1.6000 

1.6500 

1.7000 

1.7500 
1.8000 

1.8500 

1.9000 

1.9500 
2.0000 
2.0500 
2.1000 

2.1500 
2.2000 

2.2500 

2.3000 

2.3500 

2.4000 

2.4500 

2.5000 

2.5500 
2.6000 

2.6500 

2.7000 

2.7500 
2.8000 

2.8500 

2.9000 

2.9500 
3.0000 
3.0500 
3.1000 

3.1500 
3.2000 

3.2500 

3.3000 


0.604498 
0.568889 
0.533601 
0.499009 
0 . 465403 
0.433003 
0.401967 
0.372404 
0.344381 
0.317931 
0.293060 
0.269752 
0.247973 
0.227679 
0.208814 
0.191315 
0.175116 
0.160147 
0.146338 
0.133618 
0.121917 
0.111168 
0.101305 
0.092263 
0.083984 
0.076410 
0.069486 
0.063163 
0.057391 
0.052127 
0.047330 
0.042959 
0.038980 
0.035359 
0.032066 
0.029072 
0.026351 
0.023880 
0.021636 
0.019599 
0.017751 
0.016074 
0.014553 
0.013174 
0.011924 
0.010792 


0.618875 

0.648210 

0.675770 

0.701581 

0.725687 

0.748142 

0.769010 

0.788363 

0.806276 

0.822827 

0.838096 

0.852159 

0.865096 

0.876981 

0.887888 

0.897886 

0.907041 

0.915418 

0.923075 

0.930070 

0.936454 

0.942277 

0.9475-5 

0.952421 

0.956824 

0.960831 

0.964476 

0.967790 

0.970802 

0.973538 

0.976022 

0.978278 

0.980325 

0.982182 

0.983866 

0.985393 

0.986778 

0.988033 

0.989170 

0.990200 

0.991133 

0.991978 

0.992743 

0.993435 

0.994062 

0.994630 
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3.3500 

0.009766 

0.995143 

3.4000 

0.008836 

0.995608 

3.4500 

0.007994 

0.996028 

3.5000 

0.007232 

0.996409 

3.5500 

0.006541 

0.996753 

3.6000 

0.005917 

0.997064 

3.6500 

0.005351 

0.997345 

3.7000 

0.004839 

0.997600 

3.7500 

0.004376 

0.997830 

3.8000 

0.003957 

0.998038 

3.8500 

0.003578 

0.998226 

3.9000 

0.003235 

0.998397 

3.9500 

0.002925 

0.998550 

4.0000 

0.002644 

0.998690 

4.0500 

0.002391 

0.998815 

4.1000 

0.002161 

0.998929 

4.1500 

0.001954 

0.999032 

4.2000 

0.001766 

0.999125 

4.2500 

0.001597 

0.999209 

4.3000 

0.001444 

0.999285 

4.3500 

0.001305 

0.999353 

4.4000 

0.001180 

0.999415 

4.4500 

0.001067 

0.999472 

4.5000 

0.000964 

0.999522 

4.5500 

0.000872 

0.999568 

4.6000 

0.000789 

0.999610 

4.6500 

0.000713 

0.999647 

4.7000 

0.000645 

0.999681 

4.7500 

0.000584 

0.999712 

4.8000 

0.000528 

0.999740 

4.8500 

0.000478 

0.999765 

4.9000 

0.000433 

0.999788 

4.9500 

0.000392 

0.999808 

5.0000 

0.000356 

0.999827 

5.0500 

0.000323 

0.999844 

5.1000 

0.000292 

0.999859 

5.1500 

0.000266 

0.999873 

5.2000 

0.000242 

0.999886 

5.2500 

0.000221 

0.999897 

5.3000 

0.000200 

0.999908 

5.3500 

0.000185 

0.999918 

5.4000 

0.000170 

0.999927 

5.4500 

0.000157 

0.999935 

5.5000 

0.000145 

0.999942 

NUMBER  CF  POLES  EVALUATED  =  89 


204 


PROGRflM  CXMPAR  RUN  IN  DOUBLE  PRECISIO* 


X 

Y  (H-FCN) 

y  (EVAL) 

0.0500 

0.0269560422 

0.0338687172 

0.1000 

0.1047604975 

0.1147735827 

0.1500 

0.2110272898 

0.2189047501 

0.2000 

0.3266679624 

0.3300759702 

0.2500 

0.4389488040 

0.4376919836 

0.3000 

0.5402072867 

0.5352071763 

0.3500 

0.6263855876 

0.6189637207 

0.4000 

0.6958594647 

0.6873218538 

0.4500 

0.7485799015 

0.7400141770 

0.5000 

0.7854665809 

0.7776709976 

0.5500 

0.8079891806 

0.8014756091 

0.6000 

0.8178843945 

0.8129177174 

0.6500 

0.8169693788 

0.8136205051 

0.7000 

0.8070228765 

0.8052225091 

0.7500 

0.7897133043 

0.7892999138 

0.8000 

0.7665590050 

0.7673183012 

0.8500 

0.7389101583 

0.7406055613 

0.9000 

0.7079449386 

0.7103397230 

0.9500 

0.6746747268 

0.6775470467 

1.0000 

0.6399547674 

0.6431069326 

1.0500 

0.6044977968 

0.6077611209 

1.1000 

0.5688889715 

0.5721253654 

1.1500 

0.5336009937 

0.5367022861 

1.2000 

0.4990087333 

0.5018945066 

1.2500 

0.4654029237 

0.4680174730 

1.3000 

0.4330027002 

0.4353115704 

1.3500 

0.4019668816 

0.4039533060 

1.4000 

0.3724039823 

0.3740654440 

1.4500 

0.3443809949 

0.3457260524 

1.5000 

0.3179310207 

0.3189764773 

1.5500 

0.2930598403 

0.2938282943 

1.6000 

0.2697515271 

0.2702693074 

1.6500 

0.2479732076 

0.2482686787 

1.7000 

0.2276790667 

0.2277812745 

1.7500 

0.2088136949 

0.2087513166 

1.8000 

0.1913148623 

0.1911154218 

1.8500 

0.1751157994 

0.1748051079 

1.9000 

0.1601470554 

0.1597488398 

1.9500 

0.1463379965 

0.1458736804 

2.0000 

0.1336179992 

0.1331066060 

2.0500 

0.1219173875 

0.1213755376 

2.1000 

0.1111681551 

0.1106101339 

2.1500 

0.1013045089 

0.1007423864 

2.2000 

0.0922632656 

0.0917070504 

RESIDUAL 

-0.0069126750 

-0.0100130852 

-0.0078774603 

-0.0034080078 

0.0012568204 

0.0050001104 

0.0074218669 

0.0085376109 

0.0085657245 

0.0077955833 

0.0065135715 

0.0049666771 

0.0033488737 

0.0018003674 

0.0004133905 

-0.0007592962 

-0.0016954030 

-0.0023947844 

-0.0028723199 

-0.0031521652 

-0.0032633241 

-0.0032363939 

-0.0031012924 

-0.0028857733 

-0.0026145493 

-0.0023088702 

-0.0019864244 

-0.0016614617 

-0.0013450575 

-0.0010454566 

-0.0007684540 

-0.0005177803 

-0.0002954711 

-0.0001022078 

0.0000623783 

0.0001994405 

0.0003106915 

0.0003982156 

0.0004643161 

0.0005113932 

0.0005418499 

0.0005580212 

0.0005621225 

0.0005562152 
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2.2500 

0.0839841259 

0.0834419409 

0.0005421850 

2.3000 

0.0764098507 

0.0758881208 

0.0005217299 

2.3500 

0.0694863578 

0.0689899994 

0.0004963584 

2.4000 

0.0631627543 

0.0626953630 

0.0004673913 

2.4500 

0.0573913195 

0.0569553491 

0.0004359704 

2.5000 

0.0521274483 

0.0517243794 

0.0004030689 

2.5500 

0.0473295640 

0.0469600610 

0.0003695030 

2.6000 

0.0429590129 

0.0426230642 

0.0003359487 

2.6500 

0.0389799364 

0.0386769844 

0.0003029520 

2.7000 

0.0353591407 

0.0350881930 

0.0002709477 

2.7500 

0.0320659524 

0.0318256831 

0.0002402693 

2.8000 

0.0290720785 

0.0288609122 

0.0002111663 

2.8500 

0.0263514586 

0.0261676457 

0.0001838129 

2.9000 

0.0238801228 

0.0237218028 

0.0001583200 

2.9500 

0.0216360554 

0.0215013070 

0.0001347484 

3.0000 

0.0195990566 

0.0194859415 

0.0001131151 

3.0500 

0.0177506130 

0.0176572119 

0.0000934011 

3.1000 

0.0160737745 

0.0159982144 

0.0000755601 

3.1500 

0.0145530391 

0.0144935125 

0.0000595266 

3.2000 

0.0131742345 

0.0131290204 

0.0000452141 

3.2500 

0.0119244243 

0.0118918938 

0.0000325305 

3.3000 

0.0107918007 

0.0107704277 

0.0000213730 

3.3500 

0.0097656028 

0.0097539623 

0.0000116405 

3.4000 

0.0088360115 

0.0088327943 

0.0000032172 

3.4500 

0.0079940810 

0.0079980955 

-0.0000040145 

3.5000 

0.0072316960 

0.0072418379 

-0.0000101419 

3.5500 

0.0065414505 

0.0065567237 

-0.0000152732 

3.6000 

0.0059166029 

0.0059361217 

-0.0000195188 

3.6500 

0.0053510603 

0.0053740084 

-0.0000229481 

3.7000 

0.0048392420 

0.0048649142 

-0.0000256722 

3.7500 

0.0043761116 

0.0044038736 

-0.0000277620 

3.8000 

0.0039570862 

0.0039863804 

-0.0000292942 

3.8500 

0.0035779964 

0.0036083458 

-0.0000303494 

3.9000 

0.0032350765 

0.0032660609 

-0.0000309844 

3.9500 

0.0029249240 

0.0029561618 

-0.0000312378 

4.0000 

0.0026444213 

0.0026755985 

-0.0000311772 

4.0500 

0.0023907019 

0.0024216056 

-0.0000309037 

4.1000 

0.0021613239 

0.0021916766 

-0.0000303527 

4.1500 

0.0019539329 

0.0019835396 

-0.0000296067 

4.2000 

0.0017664002 

0.0017951361 

-0.0000287359 

4.2500 

0.0015968516 

0.0016246006 

-0.0000277490 

4.3000 

0.0014436486 

0.0014702432 

-0.0000265946 

4.3500 

0.0013051134 

0.0013305330 

-0.0000254196 

4.4000 

0.0011798299 

0.0012040833 

-0.0000242534 

4.4500 

0.0010667957 

0.0010896381 

-0.0000228424 

4.5000 

0.0009643209 

0.0009860600 

-0.0000217391 
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4.5500 

0.0008721016 

0.0008923190 

-0.0000202174 

4.6000 

0.0007885237 

0.0008074822 

-0.0000189585 

4.6500 

0.0007132715 

0.0007307051 

-0.0000174336 

4.7000 

0.0006450051 

0.0006612232 

-0.0000162181 

4.7500 

0.0005836256 

0.0005983440 

-0.0000147184 

4.8000 

0.0005283933 

0.0005414409 

-0.0000130476 

4.8500 

0.0004781645 

0.0004899465 

-0.0000117820 

4.9000 

0.0004330106 

0.0004433471 

-0.0000103365 

4.9500 

0.0003921517 

0.0004011779 

-0.0000090262 

5.0000 

0.0003557943 

0.0003630181 

-0.0000072238 

5.0500 

0.0003226795 

0.0003284866 

-0.0000058071 

5.1000 

0.0002919937 

0.0002972388 

-0.0000052451 

5.1500 

0.0002664966 

0.0002689625 

-0.0000024659 

5.2000 

0.0002422554 

0.0002433754 

-0.0000011200 

5.2500 

0.0002211998 

0.0002202218 

0.0000009780 

5.3000 

0.0002000859 

0.0001992705 

0.0000008154 

5.3500 

0.0001845803 

0.0001803119 

0.0000042684 

5.4000 

0.0001702609 

0.0001631567 

0.0000071042 

5.4500 

0.0001565074 

0.0001476334 

0.0000088740 

5.5000 

0.0001445091 

0.0001335868 

0.0000109223 

ESTIMATED  SSE  IS  0.0006853443 

ESTIMATED  MSE  IS  0.0000062304 

THE  MAXIMUM  ABSOLUTE  DIFFERENCE  BETWEEN  THE 
H-EUNCTION  AND  THE  ACTUAL  FUNCTION  IS 

0.0100130852 

THE  INTEGRATED  ABSOLUTE  DENSITY  DIFFERENCE 
(AN  ESTIMATE  OF  THE  AREA  BETWEEN  THE 
H-FUNCTION  AND  THE  ACTUAL  DENSITY)  IS 

0.0065133019 

EXaMPLE  3  -  ^  OF  TWO  INraTEMraWT  STftNDftRD 
UNIFORM  VRRIATES 


PROGRAM  SUMVAR  RUN  IN  DOUBLE  PRECISION 
PC®  THE  SUM  OF  2  INDEPEHDENT  H-FUNCTION  VARIATES 
INPUT  PARAMETERS  FOR  VARIABLE  NUMBEK  1  OF  THE  SUM: 
M,  N,  P,  Q  =  1011 

K,  C  =  1.0000000  1.0000000 

a,  A  =  1.00000  1.00000 

b,  B  -  0.00000  1.00000 

THE  MCMENTS  ABOUT  THE  ORIGIN  ARE: 


THE  ZEROTH  MOMENT  IS 
THE  FIRST  MOMENT  IS 
THE  SECOND  MOMENT  IS 
THE  THIRD  MOMENT  IS 


1.0000000 

0.5000000 

0.3333333 

0.2500000 
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TOE  FOURTH 
THE  FIPIH 
THE  NEXT 
THE  NEXT 
THE  NEXT 
THE  NEXT 
THE  NEXT 


MCMEirr  IS 
MQMEHT  IS 
MCMENT  IS 
MOMEHT  IS 
MQMEirr  IS 
MGMEHT  IS 
MQMEin’  IS 


0.2000000 

0.1666667 

0.1428571 

0.1250000 

0.1111111 

0.1000000 

0.0909091 


INPUT  PARAMETES^  FCX%  VARIABLE  NUMBER 
M,  N,  P,  Q  =  1011 

K,  C  =  1.0000000  1.0000000 

a,  A  =  1.00000  1.00000 

b,  B  =  0.00000  1.00000 

THE  MOMENTS  ABOUT  THE  ORIGIN  ARE: 

THE  ZEROTH  MOMENT  IS  1 

THE  FIRST  MCMEMT  IS  0 

THE  SECOND  MOMENT  IS  0 

THE  THIRD  MOMENT  IS  0 

THE  FOURTH  MOMENT  IS  0 

THE  FIFTH  MOMENT  IS  0 


2  OF  THE  SUM: 


NEXT 

NEXT 

NEXT 

NEXT 

NEXT 


MOMENT  IS 
MOMENT  IS 
MOMENT  IS 
MOMENT  IS 
MOMENT  IS 


..0000000 

>.5000000 

>.3333333 

>.2500000 

>.2000000 

>.1666667 

>.1428571 

>.1250000 

>.1111111 

>.1000000 

>.0909091 


FCR  THE  RANDOM  VARIABLE  GIVEN  AS  THE  SUM 
OF  THE  ABOVE  INDEPENDENT  H-FUNCTION 
VARIATES,  THE  MOMENTS  ABOUT  THE  CRIGIN  ARE: 


THE  ZENOTO  MOMENT  IS 
THE  FIRST  MOMENT  IS 
THE  SECOND  MOMENT  IS 
THE  THIRD  MOMENT  IS 
THE  FOURTH  MOMENT  IS 
THE  FIETO  MOMENT  IS 
THE  NBCT  MOMENT  IS 
THE  NEXT  MOMENT  IS 
THE  NEXT  MOMENT  IS 
THE  NEXT  MOMENT  IS 
THE  NEXT  MOMENT  IS 


1.0000000 

1.0000000 

1.1666667 

1.5000000 

2.0666667 

3.0000000 

4.5357143 

7.0833333 

11.3555556 

18.6000000 

31.0151515 


PRO(RAM  ILFIT  RUN  WITH  DOUBLE  PRECISION 
INPUT  WAS  OF  TYPEO 

THE  ZENOTO  MOMENT  WAS  USED  IN  THE  FIT 
RESULTS  CF  ANALYTIC  SOLUTION  - 
NUMERATOR: 
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SMRLLB(1)=  1.499999999999999667 

BIG8( 1 ) =  1 . 000000000000000000 

DElfGMINATOR: 


SMALIA( 1 )=  3 . 999999999999999334 

BIGA( 1 ) =  1 . 000000000000000000 

VALUES  CP  K  &  C  ARE: 

K=  9.027033336853110823 

C=  0.500000000000000000 


ESTIMATE)  H-FUNCTION  PARAMETEItS 
FOR  TEffi  SUM  OF  THE  INDEPEJIDEirr 


H-EUNCTION  VARIATES  ARE: 


M, 

N,  P,  Q  = 

1  0 

1 

1 

K, 

C 

9. 

0270333 

0 

.5000000 

a. 

A 

4.00000 

1.00000 

b. 

B 

1.50000 

1.00000 

THE  MQMERTS  ABOUT  THE  CPIGIN  OF  THIS 

H-EUNCTriCP  ARE: 

THE 

ZEROTH 

MCMERT 

IS 

1.0000000 

THE 

FIRST 

MOMENT 

IS 

1.0000000 

THE 

SEOCPD 

MOMENT 

IS 

1.1666667 

THE 

THIRD 

MCMERT 

IS 

1.5000000 

THE 

FOURTH 

MOMENT 

IS 

2.0625000 

THE 

FIFTH 

MOMENT 

IS 

2.9791667 

THE 

NEXT 

MCMERT 

IS 

4.4687500 

THE 

NEXT 

MOMENT 

IS 

6.9062500 

THE 

NEXT 

MCMERT 

IS 

10.9348958 

THE 

NEXT 

MCMERT 

IS 

17.6640625 

THE 

NEXT 

MCMERT 

IS 

29.0195312 

DETERMINE  P.D.F.(Z)  AND  C,D.P.(Z) 

FOR  VALUES  OF  Z  ERCM  0.0500  TO  2.0000 

WITH  STEP  SIZE  0.0500 

FOR  THE  SUM  CP  1  TERMS,  WHERE 

THE  MAXIMUM  NUMBER  CP  POLES  TO  BE  EVALUATED  IS  100 . 

CRUMP  PARAMETERS:  NUMBER  CP  CXMPLEX  VALUES  =  1001, 

PERCHTT  CP  HIGHEST  Z  VALUE  =  1.00,  AXIS  POINT  A  =  4.7785 

FORM  FCP  OVERALL  PRCPLEM  (WHERE  YJ  =  XJ**PJ): 

Z  =  Yl 

VARIATE  X  1  IS  TYPE  NUMBER  4 


INPUT  PARAMETERS  ARE  THETA  =  9.02703, 

PHI  = 

0.50000 

AND  POWER  =  1.00000 

THE  P.D.F.  FOR  VARIATE  X  1  IS  GIVER 

BY: 

1  0 

( 

4.000, 

1.000) 

9.02703  H  (  0.50000  X): 

1  1 

( 

1.500, 

1.000) 
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HIE  P.D.F.  PC«  TERM  1  OF  TOE  SUM  IS  GIVEN  BY: 
1  0 

9.02703  H  (  0.50000  Z),  WHERE 

1  1 

(BA(I),GBA(I)):  (  1.500,  1.000)  ( 
(a>(I),GCD(I)):  (  4.000,  1.000)  ( 
OOMVEKGENCE  TYPE  =  7 


D  = 

0.00  E  = 

0.00  L 

Z 

PDF(Z) 

CDFiZ) 

0.0500 

0.025842 

0.000523 

0.1000 

0.070299 

0.002876 

0.1500 

0.124083 

0.007706 

0.2000 

0.183346 

0.015375 

0.2500 

0.245633 

0.026091 

0.3000 

0.309153 

0.039958 

0.3500 

0.372518 

0.057003 

0.4000 

0.434599 

0.077189 

0.4500 

0.494465 

0.100426 

0.5000 

0.551329 

0.126585 

0.5500 

0.604526 

0.155498 

0.6000 

0.653488 

0.186967 

0.6500 

0.697733 

0.220768 

0.7000 

0.736853 

0.256655 

0.7500 

0.770506 

0.294362 

0.8000 

0.798410 

0.333610 

0.8500 

0.820340 

0.374104 

0.9000 

0.836126 

0.415541 

0.9500 

0.845645 

0.457612 

1.0000 

0.848826 

0.500000 

1.0500 

0.845645 

0.542388 

1.1000 

0.836126 

0.584459 

1.1500 

0.820340 

0.625896 

1.2000 

0.798410 

0.666390 

1.2500 

0.770506 

0.705638 

1.3000 

0.736853 

0.743345 

1.3500 

0.697733 

0.779232 

1.4000 

0.653488 

0.813033 

1.4500 

0.604526 

0.844502 

1.5000 

0.551329 

0.873415 

1.5500 

0.494465 

0.899574 

1.6000 

0.434599 

0.922811 

1.6500 

0.372518 

0.942997 

1.7000 

0.309153 

0.960042 

1.7500 

0.245633 

0.973909 

1.8000 

0.183346 

0.984625 

1.8500 

0.124083 

0.992294 

=  1.0000 
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1.9000 

0.070297  0. 

997124 

1.9500 

0.025793  0. 

999477 

2.0000 

-0.001952  0. 

999977 

NUMBEK  OF 

POLES  EVALUATED  = 

:  100 

PROGRAM  OCMPAR  RUN  IN  DOUBLE  PRECISION 

X 

Y  (H-PCN) 

y  (EVAL) 

RESIDUAL 

0.0500 

0.0258419798 

0.0500000000 

-0.0241580202 

0.1000 

0.0702990184 

0.1000000000 

-0.0297009816 

0.1500 

0.1240833016 

0.1500000000 

-0.0259166984 

0.2000 

0.1833464944 

0.2000000000 

-0.0166535056 

0.2500 

0.2456325663 

0.2500000000 

-0.0043674337 

0.3000 

0.3091534688 

0.3000000000 

0.0091534688 

0.3500 

0.3725176396 

0.3500000000 

0.0225176396 

0.4000 

0.4345990979 

0.4000000000 

0.0345990979 

0.4500 

0.4944645689 

0.4500000000 

0.0444645689 

0.5000 

0.5513288954 

0.5000000000 

0.0513288954 

0.5500 

0.6045258790 

0.5500000000 

0.0545258790 

0.6000 

0.6534882583 

0.6000000000 

0.0534882583 

0.6500 

0.6977334553 

0.6500000000 

0.0477334553 

0.7000 

0.7368531561 

0.7000000000 

0.0368531561 

0.7500 

0.7705055551 

0.7500000000 

0.0205055551 

0.8000 

0.7984095245 

0.8000000000 

-0.0015904755 

0.8500 

0.8203402272 

0.8500000000 

-0.0296597728 

0.9000 

0.8361258520 

0.9000000000 

-0.0638741480 

0.9500 

0.8456452546 

0.9500000000 

-0.1043547454 

1.0000 

0.8488263632 

1.0000000000 

-0.1511736368 

1.0500 

0.8456452546 

0.9500000000 

-0.1043547454 

1.1000 

0.8361258520 

0.9000000000 

-0.0638741480 

1.1500 

0.8203402272 

0.8500000000 

-0.0296597728 

1.2000 

0.7984095245 

0.8000000000 

-0.0015904755 

1.2500 

0.7705055551 

0.7500000000 

0.0205055551 

1.3000 

0.7368531561 

0.7000000000 

0.0368531561 

1.3500 

0.6977334553 

0.6500000000 

0.0477334553 

1.4000 

0.6534882583 

0.6000000000 

0.0534882583 

1.4500 

0.6045258790 

0.5500000000 

0.0545258790 

1.5000 

0.5513288954 

0.5000000000 

0.0513288954 

1.5500 

0.4944645689 

0.4500000000 

0.0444645689 

1.6000 

0.4345990979 

0.4000000000 

0.0345990979 

1.6500 

0.3725176396 

0.3500000000 

0.0225176396 

1.7000 

0.3091534687 

0.3000000000 

0.0091534687 

1.7500 

0.2456325660 

0.2500000000 

-0.0043674340 

1.8000 

0.1833464889 

0.2000000000 

-0.0166535111 

1.8500 

0.1240831901 

0.1500000000 

-0.0259168099 

1.9000 

0.0702967579 

0.1000000000 

-0.0297032421 

1.9500 

0.0257925005 

0.0500000000 

-0.0242074995 
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2.0000  -0.0019520985  0.0000000000  -0.0019520985 

ESTIMATED  SSE  IS  0.0920133337 

ESTIMATED  MSE  IS  0.0023003333 

'RIE  MAXIMUM  ABSOLUTE  DIFFE3^CE  BETWEEN  THE 
H-FUNCTION  AND  THE  ACTUAL  FUNCTION  IS  0.1511736368 

THE  INTEGRATED  ABSOLUTE  I^SITT  DIFTTRFNCE 
(AN  ESTIMATE  OF  THE  AREA  BETWEEN  THE 

H-PUNCTION  AND  THE  ACTUAL  DENSITY)  IS  0.0752034552 


EXAMPLE  4  -  SM  ^  3WQ  INDEPEMDEMT.  IDENTICALLY 
DISTRIBUTED  BETTA  VARIATES 


PROGRAM  SUMVAR  RUN  IN  DOUBLE  PRECISION 
PCM?  THE  SUM  CF  2  INIffiPENDENT  H-FUNCTION  VARIATES 
INPUT  PARAMET'ERS  FOR  VARIABLE  NUMBER  1  OF  THE  SUM: 
M,  N,  P,  Q  =  1011 

K,  C  =  2.0000000  1.0000000 

a,  A  =  2.00000  1.00000 

b,  B  =  0.00000  1.00000 

THE  MOMENTS  ABOUT  THE  ORIGIN  ARE: 

THE  ZEROTH  MOMENT  IS  1.0000000 

THE  FIRST  MOMENT  IS  0.3333333 

THE  SECOND  MOMENT  IS  0.1666667 

THE  THIRD  MOMENT  IS  0.1000000 

THE  FOURTH  MOMENT  IS  0.0666667 

THE  FIFTH  MOMENT  IS  0.0476190 

THE  NBCT  MOMENT  IS  0.0357143 

THE  NEXT  MOMENT  IS  0.0277778 

THE  NEXT  MOMENT  IS  0.0222222 

THE  NEXT  MOMENT  IS  0.0181818 

THE  NEXT  MOMENT  IS  0.0151515 

INPUT  PARAMETERS  FCR  VARIABLE  NUMBER  2  OF  THE  SUM: 
M,  N,  P,  Q  =  1011 

K,  C  =  2.0000000  1.0000000 

a,  A  =  2.00000  1.00000 

b,  B  =  0.00000  1.00000 

THE  MOMENTS  ABOUT  THE  CXIIGIN  ARE: 

THE  ZEROTH  MOMENT  IS  1.0000000 

THE  FIRST  MOMENT  IS  0.3333333 

THE  SECOND  MOMENT  IS  0.1666667 

THE  THIRD  MOMENT  IS  0.1000000 

THE  FOURTH  MOMENT  IS  0.0666667 
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THE  FIFTH 

MOMENT 

IS 

0.0476190 

THE 

NEXT 

MOMENT 

IS 

0.0357143 

THE 

NEXT 

MOMENT 

IS 

0.0277778 

THE 

NEXT 

MOMENT 

IS 

0.0222222 

THE 

NEXT 

MOMENT 

IS 

0.0181818 

THE 

NEXT 

MOMENT 

IS 

0.0151515 

FOR  THE  RANDOM  VARIABLE  GIVQ}  AS  THE  S(14 
OF  THE  ABOVE  INDEPEWDBNT  H-FONCTIOI 
VARIATES,  IHE  MOMENTS  ABOUT  THE  ORIGIN  ARE: 


THE 

ZEROTH  MOMENT 

IS 

1.0000000 

THE  FIRST 

MOMENT 

IS 

0.6666667 

THE 

SECOND  MOMENT 

IS 

0.5555556 

THE  THIRD 

MOMENT 

IS 

0.5333333 

THE  FOURTH  MOMENT 

IS 

0.5666667 

THE  FIFTH 

MOMENT 

IS 

0.6507937 

THE 

NESTT 

MOMENT 

IS 

0.7952381 

THE 

NEXT 

MOMENT 

IS 

1.0222222 

THE 

NEXT 

MOMENT 

IS 

1.3703704 

THE 

NEXT 

MOMENT 

IS 

1.9030303 

THE 

NEXT 

MOMENT 

IS 

2.7229437 

PROGRAM  H_FIT  RUN  WITH  DOUBLE  PRECISION 
INPUT  WAS  OF  TYPED 

THE  ZEROTH  MOMENT  WAS  USED  IN  THE  FIT 
RESULTS  OF  ANALYTIC  SOLUTIO*  - 
NUMERATOR: 

SMALLB(1)=  1.333333333333333204 

BIGB(1)=  1.000000000000000000 

DENOMINATOR: 

a<ALLA( 1 ) =  5 . 999999999999999556 

BIGA( 1 ) =  1 . 000000000000000000 

VALUES  OF  K  &  C  ARE: 

K=  302.358560865130122863 
C=  0.500000000000000000 

ESTIMATED  H-FUNCTIGN  PARAMETERS 
FOR  THE  SUM  OF  THE  INDEPENDENT 
H-FUNCTION  VARIATES  ARE: 


M, 

N,  P,  Q  =  10 

1  1 

K, 

C  =  302.3585609 

0.5000000 

a. 

A  =  6.00000 

1.00000 

b. 

B  =  1.33333 

1.00000 

THE  MOMENTS  ABOUT  THE  ORIGIN  OF  THIS 
H-FUNCTION  ARE: 

THE  ZEROTH  MOMENT  IS 


1.0000000 
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THE 

FIRST 

MOMENT 

IS 

0.6666667 

THE 

SECOND 

MOMENT 

IS 

0.5555556 

THE 

THIRD 

MOMENT 

IS 

0.5349794 

THE 

FOURTH 

MOMENT 

IS 

0.5706447 

THE 

FIFTH 

MOMENT 

IS 

0.6571060 

THE 

NEXT 

MOMENT 

IS 

0.8031296 

THE 

NEXT 

MOMENT 

IS 

1.0296533 

THE 

NEXT 

MOMENT 

IS 

1.3728711 

THE 

NEXT 

MOMENT 

IS 

1.8915113 

THE 

NEXT 

MOMENT 

IS 

2.6796410 

DETERMINE  P.D.P.(Z)  AND  C.D.F.(Z) 

FOR  VALUES  OF  Z  FRGM  0.0500  TO  2.0000 

WITH  STEP  SIZE  0.0500 

FOR  THE  SUM  OF  1  TERMS,  WHERE 

THE  MAXIMUM  NUMBER  OF  POLES  TO  BE  EVALUATED  IS  100 . 

CRUMP  PARAMETERS:  NUMBER  OF  CCMPLEX  VALUES  =  1001, 

PERCEWT  OF  HIGHEST  Z  VALUE  =  1.00,  AXIS  POINT  A  =  4.7785 

FORM  FC»  OVERALL  PROBLEM  (WHERE  YJ  =  XJ**PJ): 

Z  =  Yl 

VARIATE  X  1  IS  TYPE  NUMBER  4 

INPUT  PARAMETERS  ARE  THETA  =  302.35856,  PHI  =  0.50000 


AND  POWER 

=  1.00000 

THE  P.D.F. 

FOR  VARIATE  X 

:  1  IS  GIVEN  BY: 

1  0 

(  6 

302.35856 

H  ( 

0.50000  X): 

1  1 

(  1 

THE  P.D.F. 

FOR  TERM  1  OF  THE  SUM  IS  GIVEN 

1  0 

302.35856 

H  ( 

0.50000  Z),  WHERE 

1  1 

(BA(I),GBA(I)):  (  1.333 

i,  1.000)  ( 

(CD(I),GCD(I)):  (  6.000 

1,  1.000)  ( 

CONVERGENCE  TYPE  =  7 

D  = 

0.00  E  = 

0.00  L  =  -4.67 

Z 

PDF(Z) 

CDF(Z) 

0.0500 

0.136922 

0.003018 

0.1000 

0.313676 

0.014247 

0.1500 

0.488431 

0.034343 

0.2000 

0.648272 

0.062838 

0.2500 

0.787247 

0.098821 

0.3000 

0.902663 

0.141169 

0.3500 

0.993688 

0.188680 

0.4000 

1.060648 

0.240137 

0.4500 

1.104632 

0.294361 

0.5000 

1.127231 

0.350243 

000, 

333, 


1.000) 

1.000) 


R  =  1.0000 
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0.5500 

1.130362 

0.406760 

0.6000 

1.116150 

0.462990 

0.6500 

1.086825 

0.518123 

0.7000 

1.044659 

0.571459 

0.7500 

0.991905 

0.622412 

0.8000 

0.930756 

0.670510 

0.8500 

0.863314 

0.715383 

0.9000 

0.791555 

0.756769 

0.9500 

0.717316 

0.794498 

1.0000 

0.642271 

0.828487 

1.0500 

0.567924 

0.858736 

1.1000 

0.495600 

0.885314 

1.1500 

0.426434 

0.908349 

1.2000 

0.361374 

0.928025 

1.2500 

0.301178 

0.944568 

1.3000 

0.246416 

0.958234 

1.3500 

0.197475 

0.969306 

1.4000 

0.154565 

0.978082 

1.4500 

0.117726 

0.984864 

1.5000 

0.086842 

0.989954 

1.5500 

0.061651 

0.993643 

1.6000 

0.041760 

0.996207 

1.6500 

0.026665 

0.997899 

1.7000 

0.015767 

0.998944 

1.7500 

0.008399 

0.999535 

1.8000 

0.003848 

0.999831 

1.8500 

0.001390 

0.999954 

1.9000 

0.000326 

0.999993 

1.9500 

0.000026 

1.000000 

2.0000 

-0.000001 

1.000000 

NUMBER  OF  POLES  EVALUATED  =  100 

PROGRAM  OGMPAR  RUN  IN  DOUBLE  PRECISION 


X 

Y  (H-PCN) 

Y  (EVAL) 

RESIDUAL 

0.0500 

0.1369215248 

0.1900833333 

-0.0531618085 

0.1000 

0.3136760423 

0.3606666667 

-0.0469906244 

0.1500 

0.4884305603 

0.5122500001 

-0.0238194398 

0.2000 

0.6482720923 

0.6453333335 

0.0029387588 

0.2500 

0.7872469795 

0.7604166670 

0.0268303125 

0.3000 

0.9026634638 

0.8580000005 

0.0446634633 

0.3500 

0.9936878985 

0.9385833342 

0.0551045643 

0.4000 

1.0606483104 

1.0026666679 

0.0579816425 

0.4500 

1.1046322658 

1.0507500018 

0.0538822640 

0.5000 

1.1272305972 

1.0833333358 

0.0438972614 

0.5500 

1.1303624979 

1.1009166700 

0.0294458279 

0.6000 

1.1161501217 

1.1040000043 

0.0121501174 
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0.6500 

1.0868253832 

1.0930833388 

-0.0062579556 

0.7000 

1.0446588587 

1.0686666735 

-0.0240078148 

0.7500 

0.9919045410 

1.0312500084 

-0.0393454674 

0.8000 

0.9307563999 

0.9813333435 

-0.0505769436 

0.8500 

0.8633140228 

0.9194166789 

-0.0561026561 

0.9000 

0.7915554397 

0.8460000145 

-0.0544445748 

0.9500 

0.7173157774 

0.7615833504 

-0.0442675730 

1.0000 

0.6422707533 

0.6666666865 

-0.0243959332 

1.0500 

0.5679242735 

0.5715833504 

-0.0036590769 

1.1000 

0.4955995779 

0.4860000145 

0.0095995634 

1.1500 

0.4264335126 

0.4094166789 

0.0170168337 

1.2000 

0.3613736037 

0.3413333435 

0.0200402602 

1.2500 

0.3011776887 

0.2812500084 

0.0199276803 

1.3000 

0.2464159201 

0.2286666735 

0.0177492466 

1.3500 

0.1974750104 

0.1830833388 

0.0143916716 

1.4000 

0.1545646276 

0.1440000043 

0.0105646233 

1.4500 

0.1177258977 

0.1109166700 

0.0068092277 

1.5000 

0.0868420048 

0.0833333358 

0.0035086690 

1.5500 

0.0616509293 

0.0607500018 

0.0009009275 

1.6000 

0.0417604104 

0.0426666679 

-0.0009062575 

1.6500 

0.0266652853 

0.0285833342 

-0.0019180489 

1.7000 

0.0157674444 

0.0180000005 

-0.0022325561 

1.7500 

0.0083987678 

0.0104166670 

-0.0020178992 

1.8000 

0.0038476138 

0.0053333335 

-0.0014857197 

1.8500 

0.0013897874 

0.0022500001 

-0.0008602127 

1.9000 

0.0003256263 

0.0006666667 

-0.0003410404 

1.9500 

0.0000264951 

0.0000833333 

-0.0000568382 

2.0000 

-0.0000011508 

0.0000000000 

-0.0000011508 

ESTIMATED 

SSE  IS 

0.0358601689 

ESTIMATED  MSE  IS 

0.0008965042 

TOE  MAXIMUM  ABSOLUTE  DIFFERENCE  BEIWEEN  TOE 

H-EUNCTIC»J  AND  THE  ACTUAL 

EUNCTIOI  IS 

0.0579816425 

TOE  INTEGRATED  ABSOLUTE  MNSITY  DIFFERENCE 

(AN  ESTIMATE  OF  TOE  AREA  BETWEEN  TOE 
H-FUNCTIC»I  AND  TOE  ACTUAL  DENSITY)  IS 

0.0442126253 

EXaMPLE  5  ^  SM  TWO  INDEPEaroHfT  VRRIATES  WITH 


WEIBULL  AND  RAYLEIGH  DISTRIBUTICWS 


PROGRAM  SUMVAR  RUN  IN  DOUBLE  PRECISIOl 
PGR  THE  SUM  OF  2  INDEPENDENT  H-FUNCTIC*!  VARIATES 
INPUT  PARAMETERS  FOR  VARIABLE  NUMBER  1  OP  THE  SUM: 
M,  N,  P,  Q  =  10  0  1 
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K,  C  =  1.3195079  1.3195079 
b,  B  =  0.80000  0.20000 
THE  MOMEWTS  ABOUT  THE  ORIGIN  ARE: 


THE 

ZEROTH 

MOMENT 

IS 

1.0000000 

THE 

FIRST 

MOMENT 

IS 

0.6958418 

THE 

SEOC»n) 

MOMENT 

IS 

0.5095992 

THE 

THIRD 

MOMENT 

IS 

0.3889251 

THE 

FOURTH 

MOMENT 

IS 

0.3072421 

THE 

FIFTH 

MOMENT 

IS 

0.2500000 

THE 

NEXT 

MOMENT 

IS 

0.2087525 

THE 

NEXT 

MOMENT 

IS 

0.1783597 

THE 

NEXT 

MOMENT 

IS 

0.1555701 

THE 

NEXT 

MOMENT 

IS 

0.1382589 

THE 

NEXT 

MOMENT 

IS 

0.1250000 

INPUT  PARAMETERS  PC»  VARIABLE  NUMBER  2  OF  THE  SUM: 
M,  N,  P,  Q  =  1001 

K,  C  =  1.7320508  1.7320508 

b,  B  =  0.50000  0.50000 

THE  MOMEWTS  ABOUT  THE  C»IGIN  ARE: 


THE 

ZEROTH 

MOMENT 

IS 

1.0000000 

THE 

FIRST 

MOMENT 

IS 

0.5116634 

THE 

SEOCM) 

MOMENT 

IS 

0.3333333 

THE 

THIRD 

MOMENT 

IS 

0.2558317 

THE 

FOURTH 

MOMENT 

IS 

0.2222222 

THE 

FIFTH 

MOMENT 

IS 

0.2131931 

THE 

NEXT 

MOMENT 

IS 

0.2222222 

THE 

NEXT 

MOMENT 

IS 

0.2487252 

THE 

NEXT 

MOMENT 

IS 

0.2962963 

THE 

NEXT 

MOMENT 

IS 

0.3730879 

THE 

NEXT 

MOMENT 

IS 

0.4938272 

PC«  THE  RANDOM  VARIABLE  GIVEN  AS  THE  SUM 
OF  THE  ABOVE  INDEPEWDEWT  H-PUNCTICW 
VARIATES,  THE  MOMEWTS  ABOUT  THE  ORIGJU  ARE: 


THE 

ZENOTH  MOMENT 

IS 

1.0000000 

THE 

FIRST 

MOMENT 

IS 

1.2075051 

THE 

SEOCXO}  MOMENT 

IS 

1.5550061 

THE 

THIRD 

MOMENT 

IS 

2.1228284 

THE 

FOURTH  MOMENT 

IS 

3.0567312 

THE 

FIFTH 

MOMENT 

IS 

4.6225066 

THE 

NEXT 

MOMENT 

IS 

7.3134235 

THE 

NEXT 

MOMENT 

IS 

12.0647384 

THE 

NEXT 

MOMENT 

IS 

20.6900000 

THE 

NEXT 

MOMENT 

IS 

36.7858797 

THE 

NEXT 

MOMENT 

IS 

67.6453356 

PROGRAM  H_FIT  RUN  WITH  DOUBLE  PRECISION 

INPUT  WAS  OF  TyPEO 

DEFAULT  INITIAL  GUESS  WAS  USED 

THE  ZERCTIH  MCMEirr  WAS  USED  IN  IHE  FIT 

RESULTS  OF  ZSPOW  - 

NUMERAT(»: 

SMALLB(1)=  5.291047428929997909 

BIGB(1)=  0.623582866349619888 

DENOMINATOR: 

VALUES  OF  K  &  C  ARE: 

K=  0 . 023697857773422697 

C=  2.459933992450636764 

ENORM=  0 . 000000000000000000 

ESTIMATED  H-EUNCTIOJ  PARAMETERS 
FOR  IHE  SUM  OF  THE  INDEPENDElfr 
H-EUNCTION  VARIATES  ARE: 

M,  N,  P,  Q  =  1001 

K,  C  =  0.0236979  2.4599340 

b,  B  =  5.29105  0.62358 

THE  MCMENTS  ABOUT  THE  CKIGIN  OF  THIS 
H-FUNCTION  ARE: 


THE 

ZEROTH 

MOMENT 

IS 

1.0000000 

THE 

FIRST 

MOMENT 

IS 

1.2075051 

THE 

SECOND 

MOMENT 

IS 

1.5550061 

THE 

THIRD 

MOMENT 

IS 

2.1228284 

THE 

FOURTH 

MOMENT 

IS 

3.0568212 

THE 

FIFTH 

MOMENT 

IS 

4.6235079 

THE 

NEJCT 

MOMENT 

IS 

7.3191453 

THE 

NEXT 

MOMENT 

IS 

12.0889706 

THE 

NEXT 

MOMENT 

IS 

20.7770136 

THE 

NEXT 

MOMENT 

IS 

37.0686930 

THE 

NEXT 

MOMENT 

IS 

68.5084603 

DETERMINE  P.D.F.(Z)  AND  C.D.F.(Z) 

FOR  VALUES  OF  Z  EHOM  0.0500  TO  2.4000 

WITH  STEP  SIZE  0.0500 

FOR  THE  SUM  OF  1  TERMS,  WHERE 

THE  MAXIMUM  NUMBER  OF  POLES  TO  BE  EV2VLUATED  IS  100. 

CRUMP  PARAMETERS:  NUMBER  OF  OCMPLEX  VALUES  =  1001, 

PERCEHT  OF  HIGHEST  Z  VALUE  =  1.00,  AXIS  POINT  A  =  3.9820 

FORM  FOR  OVERALL  PROBLEM  (WHERE  YJ  =  XJ**PJ) : 

Z  =  Y1 

VARIATE  X  1  IS  TYPE  NUMBER  4 

INPUT  PARAMETERS  ARE  THETA  =  0.02370,  PHI  =  2.45993 


AND  PCWER  =  1.00000 

THE  P.D.F.  FC«  VARIATE  X  1  IS  GIVEN  BY: 

1  0 

0.02370  H  (  2.45993  X): 

0  1  (  5.291, 

TOE  P.D.F.  FC»  TERM  1  OF  TOE  SUM  IS  GIVEN  BY: 

1  0 

0.02370  H  (  2.45993  Z),  WHERE 

0  1 

(BA(I),GBA(I)):  (  5.291,  0.624)  ( 


CCXWERGENCE  TZPE  =  1 


D  = 

0.62  E  : 

Z 

PDF(Z) 

0.0500 

0.000000 

0.1000 

0.000000 

0.1500 

0.000007 

0.2000 

0.000067 

0.2500 

0.000388 

0.3000 

0.001561 

0.3500 

0.004860 

0.4000 

0.012505 

0.4500 

0.027742 

0.5000 

0.054613 

0.5500 

0.097399 

0.6000 

0.159805 

0.6500 

0.244059 

0.7000 

0.350144 

0.7500 

0.475360 

0.8000 

0.614318 

0.8500 

0.759413 

0.9000 

0.901657 

0.9500 

1.031744 

1.0000 

1.141135 

1.0500 

1.223000 

1.1000 

1.272868 

1.1500 

1.288939 

1.2000 

1.272029 

1.2500 

1.225228 

1.3000 

1.153351 

1.3500 

1.062279 

1.4000 

0.958308 

1.4500 

0.847563 

1.5000 

0.735552 

1.5500 

0.626854 

1.6000 

0.524979 

1.6500 

0.432336 

1.62  L  =  4.79  R 

CDF{Z) 

0.000000 

0.000000 

0.000000 

0.000001 

0.000011 

0.000054 

0.000202 

0.000612 

0.001579 

0.003580 

0.007305 

0.013648 

0.023651 

0.038419 

0.058986 

0.086185 

0.120521 

0.162079 

0.210484 

0.264908 

0.324137 

0.386673 

0.450859 

0.515016 

0.577563 

0.637120 

0.692578 

0.743133 

0.788296 

0.827869 

0.861907 

0.890669 

0.914559 


0.624) 


=  1.3425 
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1.7000 

0.350321 

0.934080 

1.7500 

0.279456 

0.949778 

1.8000 

0.219576 

0.962209 

1.8500 

0.170015 

0.971907 

1.9000 

0.129779 

0.979366 

1.9500 

0.097705 

0.985021 

2.0000 

0.072574 

0.989252 

2.0500 

0.053205 

0.992375 

2.1000 

0.038504 

0.994650 

2.1500 

0.027513 

0.996287 

2.2000 

0.019386 

0.997448 

2.2500 

0.013441 

0.998262 

2.3000 

0.009088 

0.998819 

2.3500 

0.005791 

0.999188 

2.4000 

0.002836 

0.999396 

NIMBES  OF  POLES  EVALUATED  =  94 


APPENDIX  B 


GRAPHICAL  DEPICTiaiS  OF  SIWS  OF  SELECTED  INDEPEWraWT 
UNIFC^,  POWER  FUNCTICW,  AND  BETA  VARIATES 


Listed  below  are  graphical  depictions  of  the  exanples  in 
Section  5.1.1  of  analytic  convolutions  performed  to  find  the 
exact  distribution  of  the  sun  of  selected  independent  random 
variables.  The  exact  distributions  of  suns  of  two  independent 
variates  with  certain  uniform,  power  function,  and  beta 
distributions  were  found  using  the  convolution  integral.  These 
graphs  nay  provide  some  insight  to  the  sum  of  independent 
variates  with  a  restricted  range.  These  sums  produce  a 
distribution  with  two  functional  forms  over  distinct  ranges  of 
the  variable  and  do  not  have  continuous  derivatives  of  all 
orders  at  x=l. 

The  two  leftmost  graphs  in  each  row  show  the  densities  of 

rightmost  graph  in  each  row  shows  the  density  of  Y=X^+X2, 
fY(y},  when  X^  and  X2  are  independent.  Each  series  of  graphs 
is  also  referenced  to  the  corresponding  equation  nxinber  from 


xj  .  The 


each  random  variable  in  the  sun. 


Section  5.1.1. 
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